Thomas Guilmeau

Emilie Chouzenoux

Víctor Elvira

ADAPTIVE SIMULATED ANNEALING THROUGH ALTERNATING R ÉNYI DIVERGENCE MINIMIZATION

Keywords: Adaptive simulated annealing, alternating Bregman proximal-gradient algorithm, Rényi divergence

Simulated annealing is a popular approach to solve nonconvex and black-box optimization problems. It consists in running a non-homogeneous Markov chain to sample from a sequence of Boltzmann probability distributions. This sequence is controlled by a cooling schedule, which governs the concentration of the mass of the Boltzmann distributions around the global minimizers. However, convergence is often slow, difficult to assess, and requires a fixed cooling schedule. We propose here a new simulated annealing algorithm with adaptive cooling schedule, which draws samples from variational approximations of the Boltzmann distributions. Our approach is theoretically sound and relies on an alternating Bregman proximal-gradient scheme minimizing a regularized Rényi divergence. Numerical experiments illustrate the performance of the method.

INTRODUCTION

In this work, we focus on the black-box minimization of a non-convex cost. In this setting, only the cost function can be evaluated, without access to any extra information (e.g. gradients). This can occur while optimizing hyper-parameters [START_REF] Watanabe | Black box optimization for automatic speech recognition[END_REF] or when dealing with categorical variables [START_REF] Abramson | Quantitative object reconstruction using Abel transform x-ray tomography and mixed variable optimization[END_REF]. Finding the global minimum in such problems requires schemes able to explore the space and avoid local minima. Many available methods are reviewed in [START_REF] Rios | Derivative-free optimization: a review of algorithms and comparison of software implementations[END_REF], including the natural evolution strategies (NES) [START_REF] Wiestra | Natural evolution strategies[END_REF] or the cross-entropy (CE) method [START_REF] Kroese | The cross-entropy method for continuous multi-extremal optimization[END_REF].

Among these methods, simulated annealing (SA) [START_REF] Haario | Simulated annealing process in general state space[END_REF] relies on Boltzmann distributions, which are probability densities constructed from the cost function and parametrized by a temperature. Low-temperature Boltzmann densities are concentrated around the global minimizers, but it is hard to sam-T.G. and E.C. acknowledge support from the ERC Starting Grant MA-JORIS ERC-2019-STG-850925. The work of V. E. is supported by the Agence Nationale de la Recherche of France under PISCES (ANR-17-CE40-0031-01), the Leverhulme Research Fellowship (RF-2021-593), and by ARL/ARO under grant W911NF-22-1-0235. ple from them. Therefore, SA schemes use decreasing temperatures to generate samples approximating these Boltzmann distributions, and thus approaching the global minimizers.

The cooling schedule (i.e., the sequence of temperatures) must be fixed in advance and decrease slowly enough to reach convergence [START_REF] Holley | Asymptotics of the spectral gap with application to the theory of simulated annealing[END_REF][START_REF] Haario | Simulated annealing process in general state space[END_REF]. In particular, the cooling schedule cannot be adapted to the quality of the samples generated in past iterations. Many adaptive SA algorithms have thus been proposed [START_REF] Molvalioglu | Meta-control of an interacting-particle algorithm for global optimization[END_REF][START_REF] Beskos | On the convergence of adaptive sequential Monte Carlo methods[END_REF] to circumvent these issues, but mostly without the theoretical guarantees of the original SA scheme [START_REF] Haario | Simulated annealing process in general state space[END_REF].

Most SA methods rely on Markov chain Monte Carlo schemes to draw samples, which allows to derive asymptotic convergence results, but can make convergence slow and difficult to assess [START_REF] Angelino | Patterns of scalable Bayesian inference[END_REF]. Sampling from parametric densities, as in the NES [START_REF] Wiestra | Natural evolution strategies[END_REF] and CE [START_REF] Kroese | The cross-entropy method for continuous multi-extremal optimization[END_REF] algorithms, provides an interesting perspective for SA schemes, as demonstrated in [START_REF] Hu | Annealing adaptive search, cross-entropy, and stochastic approximation in global optimization[END_REF].

Our recent APSA algorithm [START_REF] Guilmeau | Proximal-based adaptive simulated annealing for global optimization[END_REF] provides a methodology to adapt the cooling schedule in SA while generating parametric approximations of the Boltzmann distributions. To do so, it reformulates the problem into the minimization of the sum of the Kullback-Leibler (KL) divergence between the Boltzmann distribution and its parametric approximation, and a regularizer term. This objective function is then minimized using an alternating Bregman proximal algorithm, relying on numerical integration to adapt the temperature.

In this work, we propose a novel SA algorithm that leverages the core ideas of [START_REF] Guilmeau | Proximal-based adaptive simulated annealing for global optimization[END_REF], with a new choice of divergence, regularizer, and optimization procedure. These novel features lead to an improved algorithm which (i) is more efficient, (ii) exploits available information about the global minimum of the objective, (iii) admits theoretical guarantees, and (iv) shows good empirical performance. Specifically, we propose to use a Rényi divergence instead of the KL divergence, and we introduce a novel and more principled regularizer. The parametric distributions are chosen in an exponential family. We propose an alternating scheme involving Bregman proximal gradient steps. This allows to use the geometry of the considered families [START_REF] Nielsen | Entropies and cross-entropies of exponential families[END_REF]. We show in two benchmark cost functions that the resulting method exhibits good performance compared to existing black-box algorithms.

The rest of the paper is organized as follows. Background is provided in Section 2. We detail our method in Section 3. Numerical experiments are presented in Section 4, before concluding in Section 5.

PRELIMINARIES

Optimization problem and notation

We consider the following optimization problem

Find x * ∈ X s.t. f (x * ) = f * = min x∈X f (x). (1)
Here above, we consider a space X ⊂ R d and a cost function f : X → R. In the following, H denotes a finitedimensional Hilbert space with scalar product ⟨•, •⟩. B(X ) denotes the Borel algebra of X and given a measure ν on (X , B(X )), we denote by P(X , ν) the set of probability measures on (X , B(X )) which admits a density with respect to ν. Given a measurable function h and p ∈ P(X , ν), we write p(h) = h(x)p(x)ν(dx). For a function h, we denote its conjugate function by h * .

Exponential and Boltzmann families

Definition 1. The exponential family with sufficient statistics Γ is the family

Q = {q θ ∈ P(X , ν), θ ∈ Θ} such that q θ (x) = exp (⟨θ, Γ(x)⟩ -A(θ)) , ∀x ∈ X , (2)
with A being the log-partition function, such that Θ = dom A ⊂ H, and which reads:

A(θ) = log exp (⟨θ, Γ(x)⟩) ν(dx) , ∀θ ∈ Θ. ( 3
)
Example 1. Gaussian distributions with mean µ and covariance Σ form an exponential family [START_REF] Barndorff-Nielsen | Information and Exponential Families in Statistical Theory[END_REF], with ν being the Lebesgue measure and parameters θ = (Σ -1 µ, -1 2 Σ -1 ) ⊤ . The Boltzmann distributions, denoted by π β for β > 0, are central in SA schemes. π β gives the probability of a state x ∈ X depending on f (x) and the inverse temperature β.

Definition 2. The family B of Boltzmann distributions associated with f is the set of densities defined for β > 0 by

π β (x) = exp (-βf (x) -B(β)) , ∀x ∈ X , ( 4
)
where B is the log-partition function of π β :

B(β) = log exp (-βf (x)) ν(dx) , ∀β > 0. (5)
Remark 1. The Boltzmann family B forms an exponential family, with sufficient statistics -f and parameter β.

The mass of the Boltzmann distributions concentrates on the global minimizers of f as β increases. Indeed, lim β→+∞ π β (f ) = f * [6, Eq. (5.13)]. Therefore, points sampled from π β with high β would be close from the global minimizer. However, sampling such points is difficult.

Rényi and KL divergences

Rényi and KL divergences are widely used as discrepancy measures between probability distributions [START_REF] Li | Rényi divergence variational inference[END_REF][START_REF] Daudel | Infinite-dimensional gradient-based descent for alpha-divergence minimisation[END_REF]. Definition 3. The Rényi divergence with parameter α ∈ (0, 1) ∪ (1, +∞) and the KL divergence between p 1 and p 2 in P(X , ν) are respectively defined by

RD α (p 1 , p 2 ) = 1 α -1 log p 1 (x) α p 2 (x) 1-α ν(dx) , KL(p 1 , p 2 ) = log p 1 (x) p 2 (x) p 1 (x)ν(dx).
If one of the above is not well-defined, its value is +∞.

The KL divergence is a limiting case of Rényi divergence [START_REF] Van Erven | Rényi divergence and Kullback-Leibler divergence[END_REF], since lim α→1, α≤1 RD α (p 1 , p 2 ) = KL(p 1 , p 2 ).

PROPOSED APPROACH

A Rényi divergence minimization problem

SA methods create a sequence of Boltzmann distributions with increasing inverse temperature and approximate them with tractable distributions. In this work, we approximate the Boltzmann distributions π β ∈ B by parametric distributions q θ from an exponential family Q. Therefore, we search for parameters β, θ such that q θ and π β are close and β is sufficiently large. To achieve this goal, we consider the minimization of the following function:

J α (β, θ) = d α (β, θ) + r(β), ∀β > 0, θ ∈ Θ, (6)
where α ∈ (0, 1), d α (β, θ) = RD α (π β , q θ ). We furthermore introduce the regularizer

r(β) = B(β) + βf ϵ , ∀β > 0, (7)
where f ϵ = f * + ϵ for any small value ϵ > 0. r is convex and differentiable on R ++ . It is lower bounded as it is equal, up to an additive constant, to KL(π βϵ , π β ) where π βϵ (f ) = f ϵ . It promotes the concentration of π β around the global minimizers of f since ∇r(β) = 0 if and only if π β (f ) = f ϵ . Solving Problem (6) yields a distribution π β concentrated around the global minimizers of f , along with its approximation q θ from which samples are available. In the following, we show how to construct in an adaptive fashion a sequence {(β k , θ k )} k∈N converging to such a solution.

Proposed algorithm

We propose to minimize J α by using an alternating Bregman proximal gradient algorithm where the Bregman divergences are induced by the log-partition functions A and B. These are equivalent to the KL divergences [START_REF] Nielsen | Entropies and cross-entropies of exponential families[END_REF]:

d A (θ, θ ′ ) = KL(q θ ′ , q θ ), ∀θ, θ ′ ∈ Θ, d B (β, β ′ ) = KL(π β ′ , π β ), ∀β, β ′ > 0.
Such choices allow to decouple the updates of β and θ, so that each update reduces to a variational inference update. Indeed, in this framework, Boltzmann distributions are approached by a variational approximation from Q, and proposals from Q are also approximated by a Boltzmann distribution. Our choices also allow the use of the underlying geometries of Q and B, which often leads to more efficient schemes [START_REF] Hoffman | Stochastic variational inference[END_REF][START_REF] Raskutti | The information geometry of mirror descent[END_REF]. Finally, the available studies of [START_REF] Bauschke | Joint minimization with alternating Bregman proximity operators[END_REF][START_REF] Gheche | Proximity operators of discrete information divergences[END_REF][START_REF] Teboulle | A simplified view of first order methods for optimization[END_REF][START_REF] Ahookhosh | Multi-block Bregman proximal alternating linearized minimization and its application to orthogonal nonnegative matrix factorization[END_REF] give our method a solid theoretical background.

We are now ready to present our approach, whose iteration k ∈ N consists in the following three steps with stepsizes η k , τ k ∈ (0, 1], and η k = 1-α α η k with α ∈ (0, 1). First, we perform a Bregman gradient descent step on d α (•, θ k-1 ) within the Bregman divergence d B with step-size η k :

β k-1 2 = ∇B * (∇B(β k-1 ) -η k ∇ β d α (β k-1 , θ k-1 )
). (8) Second, we perform a Bregman proximal step on r within the divergence d B and with step-size η k :

β k = arg min β>0 r(β) + 1 η k d B (β, β k-1 2 ). (9)
Third, we perform a Bregman gradient descent step on d α (β k , •) within the divergence d A and with step-size τ k :

θ k = ∇A * (∇A(θ k-1 ) -τ k ∇ θ d α (β k , θ k-1 )). (10)
We can compute for any β > 0, θ ∈ Θ,

∇ β d α (β, θ) = α 1 -α (π β (-f ) -π (α) β,θ (-f )), (11)
∇ θ d α (β, θ) = q θ (Γ) -π (α) β,θ (Γ), (12)
where π (α)

β,θ ∈ P(X , ν) and has a density defined by

π (α) β,θ (x) ∝ exp(-f (x)) β q θ (x) α q θ (x), ∀x ∈ X . (13)
Then, using results from [START_REF] Nielsen | Entropies and cross-entropies of exponential families[END_REF][START_REF] Teboulle | A simplified view of first order methods for optimization[END_REF][START_REF] Barndorff-Nielsen | Information and Exponential Families in Statistical Theory[END_REF], we can show that Eq. ( 8)- [START_REF] Angelino | Patterns of scalable Bayesian inference[END_REF] are equivalent for every k ∈ N to

π β k-1 2 (f ) = (1 -η k ) π β k-1 (f ) + η k π (α) β k-1 ,θ k-1 (f ), (14) π β k (f ) = 1 1 + η k π β k-1 2 (f ) + η k 1 + η k f ϵ , (15)
q θ k (Γ) = (1 -τ k )q θ k-1 (Γ) + τ k π (α) β k ,θ k-1 (Γ). (16)
Let us now state our convergence result for the proposed iterative scheme, using denominations from [START_REF] Barndorff-Nielsen | Information and Exponential Families in Statistical Theory[END_REF].

Proposition 1. Let ϵ > 0. Assume that B and Q are regular and full, that τ k , η k ∈ (0, 1] for all k ∈ N, and that the updates are well-defined. Then, the sequence

{J α (β k , θ k )} k∈N is decreasing. If further, B ⊂ Q, then a fixed point of the algorithm (β * , θ * ) is such that π β * (f ) = q θ * (f ) = f ϵ .
Proposition 1 establishes a monotonic decrease of the objective J α . The second part of the result applies in particular when f is quadratic and Q is the Gaussian family.

Our Rényi-based adaptive SA (RASA) scheme

We now discuss the practical implementation of our algorithm for black-box optimization, meaning that only evaluations of f are available. We propose to approximate π (α) β,θ and π β using adaptive importance sampling [START_REF] Bugallo | Adaptive importance sampling: The past, the present, and the future[END_REF]. To this end, consider θ ∈ Θ and denote by X θ a set of N samples from q θ . We define for β > 0, α > 0 the empirical distributions

π (β,α) X θ ∝ x∈X θ exp(-f (x)) β q θ (x) α δ x . ( 17
)
The distribution above approximates π

(α) β,θ , while π (β,1) X θ
approximates π β . This yields Algorithm 1, which approximates the idealized iterative algorithm described by Eq. ( 14)-( 16), with a precision depending on the sample size N .

Algorithm 1: RASA algorithm input: Choose β 0 > 0, θ 0 ∈ Θ, N ∈ N, and {τ k , η k } k∈N in (0, 1]. Set η k = 1-α α η k for every k ∈ N. for k = 1, 2, . . . do 1 Sample N points X θ k-1 = {x n } N n=1 from q θ k-1 . 2 Set f k = min{f (x), x ∈ X θ l , 0 ≤ l ≤ k -1}. 3 If k = 1, set π (β k-1 ,1) X θ k-2 (f ) = π (β k-1 ,1) X θ k-1 (f ).
4 Compute the following empirical moment:

π β k-1 2 (f ) = (1 -η k ) π (β k-1 ,1) X θ k-2 (f )+η k π (β k-1 ,α) X θ k-1 (f ). 5 Set β k such that π (β k ,1) X θ k-1 (f ) = 1 1 + η k π β k-1 2 (f ) + η k 1 + η k f k .
6 Adapt the proposal by

q θ k (Γ) = (1 -τ k )q θ k-1 (Γ) + τ k π (β k ,α) X θ k-1 (Γ). Note that β -→ π (β,1) X θ k-1
(f ) is a decreasing continuous function. Thus, β k can be found by a dichotomy search.

Discussion

In [START_REF] Guilmeau | Proximal-based adaptive simulated annealing for global optimization[END_REF], KL divergences and quadratic regularizers are used. In contrast, our scheme relies on Rényi divergences benefiting from a useful symmetry property [17, Proposition 2] and a novel regularizer, leading to simpler β updates.

The adaptation step of θ of the MARS algorithm [START_REF] Hu | Annealing adaptive search, cross-entropy, and stochastic approximation in global optimization[END_REF] is recovered in Algorithm 1 with α = 1. However, a fixed cooling schedule was used in [START_REF] Hu | Annealing adaptive search, cross-entropy, and stochastic approximation in global optimization[END_REF], while ours is fully adaptive.

Non-linear importance weights arise in Algorithm 1 from the densities π (α) β,θ and the approximations π (β,α) X θ . Such weights are known to prevent weight degeneracy issues [START_REF] Koblents | A population Monte Carlo scheme with transformed weights and its application to stochastic kinetic models[END_REF].

NUMERICAL SIMULATIONS

Test problems and compared algorithms

We use two benchmarks with minimum value f * reached at a unique x * . Both are defined for any x ∈ R d , where d ∈ N. We consider a banana-shaped Rosenbrock function and a multimodal Rastrigin function, defined respectively by

f (x -x * ) = d-1 i=1 10(x i+1 -1 -(x i -1) 2 ) 2 + x 2 i + f * , f (x -x * ) = 4d + d i=1 0.4x 2 i -4 cos(2πx i ) + f * .
The minimum values f * are sampled uniformly in [-1, 1], and the two functions are translated so that they reach the value f * at x * uniformly sampled in [-1, 1] d . We compare our method with other schemes updating approximating densities using sampling. Namely, we consider two SA schemes, the APSA scheme [START_REF] Guilmeau | Proximal-based adaptive simulated annealing for global optimization[END_REF] and the MARS algorithm [START_REF] Hu | Annealing adaptive search, cross-entropy, and stochastic approximation in global optimization[END_REF], as well as the cross-entropy (CE) algorithm [START_REF] Kroese | The cross-entropy method for continuous multi-extremal optimization[END_REF].

All the algorithms are run with step-size τ k = 0.5 k+1 for the parameters updates, and step-size η k = 0.9 for the eventual temperature updates. The CE algorithm is run using the proportion ρ = 0.5. We simulated N = 100 samples at each iteration k ∈ N from Gaussian proposals with means µ k . The algorithms are initialized with mean µ 0 uniformly sampled in [-5, 5] d , covariance Σ 0 = 10I, and initial temperature β 0 = 0.1. At each iteration, the dichotomy search in RASA is performed in the interval [0.1β k-1 , 1.5β k-1 ].

Assessing our adaptive cooling schedule

We first compare the cooling schedules of the APSA, MARS, and RASA algorithms in dimension d = 2. Indeed, the APSA scheme is implemented with numerical integration steps which become intractable in higher dimensions. Figure 1 shows that the three methods have very different cooling schedules. As observed in [START_REF] Guilmeau | Proximal-based adaptive simulated annealing for global optimization[END_REF], the cooling schedule of the APSA algorithm reaches a high stationary value very fast. The MARS algorithm follows a logarithmic cooling schedule, which is the slowest of the three. In contrast, the cooling schedule of the APSA algorithm increases in a convex fashion, much faster than the logarithmic cooling schedule and without the stagnation of the APSA schedule.

Optimization performance in high dimension

We now compare our RASA algorithm with the MARS and CE algorithms for d = 50, where the task is to minimize the benchmark objective f . Note that the APSA method is intractable in this high-dimensional setting. We can see in Fig. 2 that RASA outperforms the MARS and CE algorithms on both benchmarks. This demonstrates the positive impact on optimization performance of our adaptive cooling schedule. Note that contrary to the CE algorithm, which only uses a ranking of the samples, the RASA and MARS algorithms use each sample's cost function value.

We also note that the lowest value of α yields the best performance. Indeed, transforming the weights as in Eq. ( 17) reduces weight degeneracy [START_REF] Koblents | A population Monte Carlo scheme with transformed weights and its application to stochastic kinetic models[END_REF], which may lower the approximation error arising from sampling with a limited number of points in high dimension. Also, low values of α lead to approximations q θ fitting the mode of π β [START_REF] Li | Rényi divergence variational inference[END_REF]. This shows the interest of the additional parameter α, coming from the use of a Rényi divergence instead of a KL divergence as in [START_REF] Hu | Annealing adaptive search, cross-entropy, and stochastic approximation in global optimization[END_REF].

CONCLUSION

In this work, we presented a novel SA algorithm for blackbox optimization. It approaches Boltzmann distributions by variational approximations and treats the Boltzmann distributions as variational approximations themselves. This creates an adaptive cooling schedule while avoiding some drawbacks of the commonly used Markov chains schemes. This procedure is interpreted as an alternating Bregman proximalgradient algorithm, allowing to gain novel theoretical insights. Numerical experiments illustrate the good behavior of the resulting method. This work also opens potential research tracks. On a theoretical level, the links between our practical implementation and its exact counterpart remain to be established precisely. The use of more complex approximating distributions, such as mixtures, can also be investigated, in particular in relation to multi-modal objectives.

Fig. 1 .

 1 Fig. 1. Plots of the quantities β k averaged over 500 runs for the Rastrigin (left) and Rosenbrock (right) functions in dimension d = 2.

Fig. 2 .

 2 Fig. 2. Plots of the quantity f (µ k ) -f * averaged over 500 runs for the Rastrigin (left) and Rosenbrock (right) functions in dimension d = 50.


Émilie Chouzenoux

Víctor Elvira

GRAPHIT: ITERATIVE REWEIGHTED 1 ALGORITHM FOR SPARSE GRAPH INFERENCE IN STATE-SPACE MODELS

Keywords: State-space model, graphical inference, sparsity, iterative reweighted scheme, MM algorithm

State-space models (SSMs) are a common tool for modeling multi-variate discrete-time signals. The linear-Gaussian (LG) SSM is widely applied as it allows for a closed-form solution at inference, if the model parameters are known. However, they are rarely available in real-world problems and must be estimated. Promoting sparsity of these parameters favours both interpretability and tractable inference. In this work, we propose GraphIT, a majorization-minimization (MM) algorithm for estimating the linear operator in the state equation of an LG-SSM under sparse prior. A versatile family of non-convex regularization potentials is proposed. The MM method relies on tools inherited from the expectationmaximization methodology and the iterated reweighted-l1 approach. In particular, we derive a suitable convex upper bound for the objective function, that we then minimize using a proximal splitting algorithm. Numerical experiments illustrate the benefits of the proposed inference technique.

INTRODUCTION

Signal processing applications often involve dealing with discrete-time multi-variate series that need to be processed to perform statistical tasks such as prediction. State-space models (SSMs) relate the observation process with a hidden state that evolves over time. This versatile approach allows to describe complex phenomena through the dynamics of the hidden state. It also allows for the incorporation physical knowledge about the system and to retrieve relevant structure that can be inferred when processing the data. There exists abundant literature in graphical models for multi-variate time series, e.g., [START_REF] Eichler | Graphical modelling of multivariate time series[END_REF][START_REF] Bach | Learning graphical models for stationary time series[END_REF][START_REF] Barber | Graphical models for time-series[END_REF][START_REF] Ioannidis | Semi-blind inference of topologies and dynamical processes over dynamic graphs[END_REF], and more recently also within SSMs [START_REF] Chouzenoux | GraphEM: EM algorithm for blind Kalman filtering under graphical sparsity constraints[END_REF][START_REF] Elvira | Graphical inference in linear-Gaussian state-space models[END_REF], including also fully probabilistic approaches [START_REF] Cox | Parameter estimation in sparse linear-Gaussian state-space models via reversible jump Markov Chain Monte Carlo[END_REF][START_REF] Cox | Sparse Bayesian estimation of parameters in lineargaussian state-space models[END_REF]. There are multiple applications of such graphical representa-V.E. acknowledges support from the Agence Nationale de la Recherche of France under PISCES (ANR-17-CE40-0031-01) project, the Leverhulme Research Fellowship (RF-2021-593), and by ARL/ARO under grant W911NF-22-1-0235. E.C. acknowledges support from the European Research Council under Starting Grant MAJORIS ERC-2019-STG-850925.

tions, e.g., in biology [START_REF] Pirayre | BRANE Clust: Clusterassisted gene regulatory network inference refinement[END_REF][START_REF] Luengo | Hierarchical algorithms for causality retrieval in atrial fibrillation intracavitary electrograms[END_REF], social network analysis [START_REF] Ravazzi | Learning influence structure in sparse social networks[END_REF], or neuroscience [START_REF] Richiardi | Machine learning with brain graphs[END_REF].

The graphical approach in SSMs benefits from the promotion of sparsity with, at least, the following three advantages. First, sparsity enhances interpretability, which is highly desired for instance in models where the hidden state has a physical interpretation. Second, it reduces the effective dimension of the parameters, alleviating the curse of dimensionality during estimation. Third, it is possible to promote other properties that are complementary with sparsity, allowing to introduce useful prior knowledge on the system, e.g., stability of the hidden process and spectral constraints [START_REF] Chiu | Low-rank constraints for fast inference in structured models[END_REF].

In [START_REF] Chouzenoux | GraphEM: EM algorithm for blind Kalman filtering under graphical sparsity constraints[END_REF][START_REF] Elvira | Graphical inference in linear-Gaussian state-space models[END_REF], we introduced GraphEM, an expectationmaximization (EM) method for parameter estimation in linear Gaussian (LG) SSMs. Specifically, we addressed the challenging task of estimating the transition matrix that encodes the Markovian dependencies in the evolution of the multi-variate state. GraphEM is based on a novel perspective that relates the transition matrix to the adjacency matrix of a directed graph, encoding causal relationship among state dimensions in the Granger-causality sense. The M-step in the EM procedure relies on a proximal splitting, benefiting from sound convergence guarantees. GraphEM thus extends the well-established EM technique for maximum likelihood estimation of LG-SSM parameters [START_REF] Roweis | A unifying review of linear Gaussian models[END_REF][START_REF] Ghahramani | Parameter estimation for linear dynamic systems[END_REF][START_REF] Shumway | An approach to time series smoothing and forecasting using the EM algorithm[END_REF][START_REF] Cappe | Inference in Hidden Markov Models[END_REF] into a maximum a posteriori estimator that suitably encompasses sophisticated prior terms. However, it is limited to convex regularizers. In the literature of sparse graph inference, non-convex penalties such as SCAD [START_REF] Fan | Variable selection via nonconcave penalized likelihood and its oracle properties[END_REF], MCP [START_REF] Zhang | Nearly unbiased variable selection under minimax concave penalty[END_REF], CEL0 [START_REF] Soubies | A continuous exact l0 penalty (cel0) for least-squares regularized problem[END_REF] have shown to be more suited than their convex surrogate 1 norm (i.e., Lasso) to properly reconstruct very sparse graphs [START_REF] Ying | Nonconvex sparse graph learning under Laplacian constrained graphical model[END_REF][START_REF] De Miranda Cardoso | Graphical models in heavytailed markets[END_REF][START_REF] Williams | Beyond Lasso: a survey of nonconvex regularization in Gaussian graphical models[END_REF]. These penalties can be viewed as continuous non-convex approximations to the quasi-norm 0 . However, non-convexity prevents the use of standard optimization tools. An efficient strategy is to resort to iterative reweighted (IR) schemes [START_REF] Wipf | Iterative reweighted 1 and 2 methods for finding sparse solutions[END_REF][START_REF] Huber | Robust Statistics[END_REF], which are based on the majorizationminimization (MM) principle [START_REF] Sun | Majorization-minimization algorithms in signal processing, communications, and machine learning[END_REF] to recast the complicated objective function into a sequence of tractable upper bounds. IR 1 algorithms have been widely applied to several signal processing problems [START_REF] Daubechies | Iteratively reweighted least squares minimization for sparse recovery[END_REF][START_REF] Chartrand | Iteratively reweighted algorithms for compressive sensing[END_REF][START_REF] Ochs | On iteratively reweighted algorithms for nonsmooth nonconvex optimization in computer vision[END_REF][START_REF] Carrillo | Sparsity Averaging Reweighted Analysis (SARA): a novel algorithm for radio-interferometric imaging[END_REF][START_REF] Sun | Graphical nonconvex optimization via an adaptive convex relaxation[END_REF][START_REF] Ying | Does the 1-norm learn a sparse graph under Laplacian constrained graphical models?[END_REF], but, to our knowledge, remain unstudied within SSMs.

In this paper, we propose a novel method called GraphIT that combines the EM technique with the powerful IR 1 strategy for the estimation of model parameters in SSMs. More precisely, GraphIT is able to provide an estimation of the LG-SSM transition matrix for a wide class of non-convex sparsity enhancing priors. Thanks to the IR 1 framework, GraphIT algorithm benefits from simple and not expensive updates. Experimental results on several synthetic datasets allow us to assess the ability of our method for improving the recovery rate for sparse transition matrices. The paper is structured as follows. In Section 2, we describe the model, the filtering/smoothing algorithms, and revisit the GraphEM framework. The novel GraphIT algorithm is presented in Section 3. The paper ends with numerical results in Section 4 and some concluding remarks in Section 5.

BACKGROUND

Inference in linear state-space model

We consider the LG-SSM described by

x k = Ax k-1 + q k , y k = Hx k + r k , (1)
for k = 1, . . . , K, where x k ∈ R Nx and y k ∈ R Ny are the hidden state and observation at time k, respectively,

A = R Nx×Nx , H = R Ny×Nx , and {q k } K k=1 ∼ N (0, Q) and {r k } K k=1 ∼ N (0, R) the i.i.d.
state and observation noise processes, respectively. The distribution of the initial state is x 0 ∼ N (x 0 ; µ 0 , Σ 0 ) with known µ 0 and Σ 0 .

In SSMs, the filtering problem consists in the computation of the so-called filtering distributions, p(x k |y 1:k ), where y 1:k {y j } k j=1 . In the LG-SSM of Eq. ( 1), the Kalman filter (KF) exactly computes the filtering distributions in a sequential and efficient manner, yielding p(x k |y 1:k ) = N (x k |µ k , Σ k ), i.e., for every k = 1, . . . , K, it computes the mean µ k and covariance Σ k of each Gaussian filtering distribution [START_REF] Kalman | A new approach to linear filtering and prediction problems[END_REF]. The smoothing problem aims at computing p(x k |y 1:K ), i.e., the posterior probability density functions (pdfs) conditioning on the whole set of observations instead, for all k = 1, . . . , K. In the LG-SSM, the smoothing pdfs are also available in a closed-form expression (again, Gaussian pdfs) by the Rauch-Tung-Striebel (RTS) smoother,

p(x k |y 1:K ) = N (x k |µ s k , Σ s k ).
Both KF and RTS smoother require to know the model parameters µ 0 , Σ 0 , A, H, Q, and R. The estimation of the model parameters in the model above has been a subject of intense research, e.g., recently [START_REF] Sharma | Recurrent dictionary learning for state-space models with an application in stock forecasting[END_REF][START_REF] Sharma | Blind kalman filtering for short-term load forecasting[END_REF] (see [START_REF] Sarkka | Bayesian Filtering and Smoothing[END_REF]Ch. 12] for a review).

GraphEM to estimate the transition matrix

GraphEM [START_REF] Chouzenoux | GraphEM: EM algorithm for blind Kalman filtering under graphical sparsity constraints[END_REF][START_REF] Elvira | Graphical inference in linear-Gaussian state-space models[END_REF] is an EM-based algorithm that we recently introduced for the maximum a posteriori (MAP) estimation of A under log-concave prior distribution p(A). GraphEM aims at finding the minimum of L(A) L 0 (A) + L 1:K (A) with L 0 (A) -log p(A) and L 1:K (A) -log p(y 1:K |A). The matrix A is interpreted in GraphEM to encode a directed graph among the state dimensions. Prior p(A) can encode structural information on such graph (e.g., sparsity). This graph plays an important role in the SSM interpretation, with links to Granger-causality models (see a discussion in [START_REF] Elvira | Graphical inference in linear-Gaussian state-space models[END_REF]). The likelihood of A associated to the model in Eq. ( 1) is

L 1:K (A) = K k=1 1 2 log|2πS k |+ 1 2 z k S -1 k z k , (2)
where z k = y k -HAµ k-1 and S k is the covariance matrix of the predictive distribution p(y k |y 1:k-1 , A) = N (y k ; HAµ k-1 , S k ), both being obtained by the KF, run for a given A (see [START_REF] Sarkka | Bayesian Filtering and Smoothing[END_REF]Section 4.3]). The direct minimization of L is difficult due to (i) the implicit form of L 1:K , and (ii) the possibly non-differentiability of L 0 . The GraphEM approach implements an EM strategy, generalizing the one in [START_REF] Shumway | An approach to time series smoothing and forecasting using the EM algorithm[END_REF], by relying on successive upper bounds of L 1:K , that leads to obtain a sequence of tractable inner problems.

Let us consider KF/RTS outputs for a given A ∈ R Nx×Nx . Denote, for every k ∈ {1, . . . , K},

G k = Σ k (A ) (A Σ k (A ) + Q) -1 , and set Ψ = K k=1 Σ s k + µ s k (µ s k ) , Φ = K k=1 Σ s k-1 + µ s k-1 (µ s k-1 ) , and ∆ = K k=1 Σ s k G k-1 + µ s k (µ s k-1
) . Then, according to [36, Chap. 12],

(∀A ∈ R Nx×Nx ) L 1:K (A) ≤ -p(x 0:K |y 1:K , A ) log p(x 0:K , y 1:K |A)dx 0:K + ct /A = 1 2 tr Q -1 (Σ -∆A -A∆ + AΦA ) + ct /A , (3)
where tr is the trace operator, and ct /A is a constant term not depending on A, such that equality holds at A = A . The GraphEM algorithm iterates alternating between an Estep, building the right term in [START_REF] Barber | Graphical models for time-series[END_REF], and an M-step minimizing the sum of this term with L 0 , using a proximal splitting scheme. In the simple case of a Laplace prior (i.e., 1 norm penalty), the minimization in the M-step can be performed by a Douglas-Rachford algorithm, leading to the simpler form of GraphEM in [START_REF] Chouzenoux | GraphEM: EM algorithm for blind Kalman filtering under graphical sparsity constraints[END_REF], while a more sophisticated primal-dual splitting technique is proposed in [START_REF] Elvira | Graphical inference in linear-Gaussian state-space models[END_REF] for dealing with a generic convex penalty. Convergence guarantees for the resulting GraphEM iteration are discussed in [START_REF] Elvira | Graphical inference in linear-Gaussian state-space models[END_REF].

PROPOSED GRAPHIT APPROACH

Considered class of penalties

Our study focuses on the class of non-convex sparsifying penalties that have been studied for instance in [START_REF] Gasso | Recovering sparse signals with a certain family of nonconvex penalties and dc programming[END_REF] (see also [START_REF] Pan | Iteratively reweighted 1-penalized robust regression[END_REF] in the context of robust estimation). We define γ exp(-λ|u|)

L 0 (A) = 1≤i,j≤Nx ρ(|A i,j |), (4)
MCP γ|u|-u 2 2λ if |u|≤ λγ λγ 2 2 otherwise , with λ > 0 γ -|u| λ if |u|≤ λγ 0 otherwise SCAD      γ|u| if |u|≤ γ -γ 2 -2aγ|u|+u 2 2(a-1) if γ < |u|≤ aγ (a+1)γ 2 2 otherwise , with a > 2      γ if |u|≤ γ γ -2|u|-(2a-1)γ 2(a-1) if γ < |u|≤ aγ 0 otherwise
Table 1: Example of potential functions and their derivatives with hyper-parameter γ > 0 such that ρ (0 + ) = γ.

where A i,j denotes the (i, j) entry of A, and ρ : R → [0, +∞) is a potential satisfying:

1. ρ is continuous, and it is increasing on [0, +∞) with ρ(0) = 0;

2. ρ is differentiable on (0, +∞), with derivative (when defined) denoted ρ ;

3. ρ is decreasing on (0, +∞) and lim u→0 + ρ (u) = γ ∈ (0, +∞).

A large class of potential functions satisfy the above requirements. Some relevant examples are listed in Table 1. Except for SCAD, in our examples, epi-convergence of u → ρ(|u|) to the discrete potential u → δ u =0 (equals 0 for a zero-valued input; 1 otherwise) is established for λ → 0 + . Then, the prior in (4) can be seen as an approximation of γ 0 . Fig. 1 shows the plots for a subset of potentials, with γ = 1. We also show 1 and 0 penalties for comparison. The considered class of priors is continuous, but still non-convex approximations of 0 .

IR 1 approach

The properties assumed on the potential ρ yield a key property, that is at the core of the so-called IR 1 optimization approach [START_REF] Zou | One-step sparse estimates in nonconcave penalized likelihood models[END_REF]. Let v ∈ R. Then,

(∀u ∈ R) ρ(|u|) ≤ ρ (|v|)(|u|-|v|) + ρ(|v|), (5)
where equality holds at u = v. A proof for this result can be found in [START_REF] Ochs | On iteratively reweighted algorithms for nonsmooth nonconvex optimization in computer vision[END_REF]. The following convex upper bound is then obtained, for (4), at some A ∈ R Nx×Nx :

L 0 (A) ≤ Ω(A ) A 1 +ct /A , (6)
where ct /A is such that equality holds at A = A . Hereabove, Ω(A ) ∈ [0, +∞) Nx×Nx is a matrix of weights defined as Ω i,j = ρ (|A i,j |) for every (i, j). Moreover, is the Hadamard (i.e., element-wise) product. We provide in Table 1 Ω(A ) • 1 +f . When f is a least square discrepancy, the inner problem in IR 1 can be solved using one of the several efficient approaches proposed in the literature for penalized least squares under Lasso penalty [START_REF] Combettes | Fixed point strategies in data science[END_REF][START_REF] Daubechies | An iterative thresholding algorithm for linear inverse problems with a sparsity constraint[END_REF][START_REF] Beck | A fast iterative shrinkage-thresholding algorithm for linear inverse problems[END_REF].

GraphIT algorithm

In the following, we describe the GraphIT algorithm, which is summarized in Algorithm 1. GraphIT is an MM algorithm, to search for the MAP estimate of A within our SSM (i.e., the minimizer of L), under the prior defined in ( 4). It has an iterative structure, which alternates over majorization and minimization steps. The difficulty lies in the construction of the majorization approximation and in the required optimization algorithm to minimize it.

Algorithm 1 GraphIT algorithm

Initial A (0) and penalty parameters. Precision ε > 0. For i = 1, 2, . . .

Majorization:

Run KF and RTS smoother using A (i-1) and construct

(∀A ∈ R Nx×Nx ) Q(A; A (i-1) ) = 1 2 tr Q -1 (Σ-∆A -A∆ +AΦA ) + Ω(A (i-1) ) A 1 . (7) Minimization:
Solve A (i) = argmin A Q(A; A (i-1) ) using DR solver.

Stopping condition:

If

A (i) -A (i-1) F ≤ ε A (i-1)
F , stop the recursion.

Given an estimation, A (i-1) , coming from the previous iteration, we construct a specific majorization function Q(A; A (i-1) ) defined in Eq. ( 7) that combines (i) an EM-based majorization involving outputs from running the KF/RTS with transition matrix A (i-1) , and (ii) a convex upper bound for our prior term in Eq. (4) (coming from IR 1 approach). Using (3) and ( 6), we can show that the function in Eq. ( 7) is such that, for every A ∈ R Nx×Nx and every iteration i ∈ N, L(A) ≤ Q(A; A (i-1) ) + ct /A (heragain, the constant term is set so as to have equality at A = A (i-1) ). In practice, the minimization of the above upper bound is done through the Douglas-Rachford (DR) algorithm [START_REF] Combettes | A Douglas-Rachford splitting approach to nonsmooth convex variational signal recovery[END_REF] initialized at the current iterate A (i-1) . As a virtue of the MM construction, the sequence {L(A (i) )} i∈N is a decreasing sequence converging to a finite limit.

EXPERIMENTAL EVALUATION

We illustrate the performance of the proposed GraphIT method on a time series simulated according to model (1), with N x = N y ∈ {8, 16}. Sparse matrices A are randomly generated, with support size (i.e., number of edges in the associated graph) S ∈ {4, 8, 16}. For each of the six scenarios, we randomly select S indexes for which we set the associated entry in A equals to a realization of an i.i.d. standard normal distribution and with the N 2

x -S other entries being set to zero. The non-zero entries of A are scaled so that the spectral norm of A is less than one to guarantee the stability of the Markov process. We finally set H = Id,

K = 10 3 and Q = σ 2 Q Id, R = σ 2 R Id, Σ 0 = σ 2 0 Id with (σ Q , σ R , σ 0 ) = (10 -1 , 10 -1 , 10 -4 ).
For each of the six datasets, we ran GraphIT using logsum penalty (Table 1, first row) and precision ξ = 10 -3 (with a maximum of 50 iterations). We also display the results of GraphEM using 1 penalty, and of the maximum likelihood estimator computed by an EM algorithm (MLEM) [START_REF] Ghahramani | Parameter estimation for linear dynamic systems[END_REF][START_REF] Shumway | An approach to time series smoothing and forecasting using the EM algorithm[END_REF]. All methods are initialized with a matrix with entries A i,j = (10 -1 ) |i-j| , rescaled to have spectral norm equals to 0.99. The methods are compared in terms of relative mean square error (RMSE) on A, accuracy and F1 scores for the graph edge detection, using a threshold of 10 -10 on the entries of matrix A. All results are averaged on 50 realizations. The regularization parameters of GraphIT and GraphEM are finetuned through a grid-search strategy to minimize the RMSE on a single realization.

The results are summarized in Table 2. MLEM does not include sparsity prior on A, which explains its poor performance, especially in terms of edge detectability. GraphEM and GraphIT both reach better scores than MLEM. It is noticeable that, the sparser the graph (i.e., the smaller S), the better results for GraphIT. The superiority of GraphIT over its competitors on all datasets is remarkable both in terms of RMSE and detection metrics. In particular, despite integrating a sparsifying 1 penalty in its estimator, GraphEM struggles for very sparse graphs, leading to an increase in RMSE and a lowered accuracy. In contrast, GraphIT performs well for all metrics. When the level of sparsity decreases (i.e., S increases), both GraphEM and GraphIT tend to perform similarly. We also include in the last column of Table 2, the averaged computing times over 50 realizations for each methods. All methods take a similar running time (code is run in Matlab 2021a on an 11th Gen Intel(R) Core(TM) i7-1185G7 3.00GHz with 32 Go RAM).

We finally display two examples of graph reconstruction in Fig. 2. We compare the true graph (left), with GraphEM (middle) and GraphIT (right) estimates in both setups, i.e., (N x , S) = (8, 16) (top) and (N x , S) = (16, 16) (bottom). We can see the superior ability of GraphIT to recover the graphs shape and weights in the case of two sparse graphs.

CONCLUSION

In this paper, we have proposed GraphIT, an MM algorithm for the estimation of a sparse transition matrix describing hidden state interactions in an LG-SSM. A family of non-convex regularization term is considered to enforce the sparsity of the sought matrix, here interpreted as related to the adjacency matrix of a directed graph. The novel method exploits majorization properties inherited from both the EM framework and the iterated reweighted 1 methodology. The resulting convex upper bounds can be efficiently minimized through a proximal splitting solver. Through numerical results on controlled scenarios, we illustrate the great ability of the method to properly retrieve sparse transition matrices.

γ 1+λ 2 u 2

 2 γλ(log(|u|+λ) -log(λ)), with λ > 0 γλ(|u|+λ) -1 atan γ λ atan(λ|u|), with λ > 0 Mangasarian γ λ (1 -exp(-λ|u|))

  the expression for ρ (|•|) for our examples. The IR 1 technique amounts at minimizing L 0 + f , where f a simple function (typically a quadratic term), by iteratively solving the tractable problem of the minimization of

Fig. 1 :

 1 Fig. 1: Example of functions ρ(|•|), with λ = 0.5 and γ = 1.

Fig. 2 :

 2 Fig. 2: True graph (left), GraphEM (middle), and GraphIT (right) estimates for (N x , S) = (8, 16) (top) and (N x , S) = (16, 16) (bottom).

Table 2 :

 2 Results and computing times for GraphIT, GraphEM, and MLEM.

	(N x , S)	method	RMSE accur.	F1	Time (s.)
		GraphIT	0.185	0.964 0.776	1.523
	(8,4)	GraphEM	0.244	0.920 0.635	1.186
		MLEM	0.401	0.063 0.118	1.862
		GraphIT	0.161	0.976 0.901	1.907
	(8,8)	GraphEM	0.214	0.936 0.784	3.467
		MLEM	0.345	0.125 0.222	1.774
		GraphIT	0.190	0.931 0.861	1.663
	(8,16)	GraphEM	0.194	0.865 0.772	3.506
		MLEM	0.248	0.250 0.400	2.040
		GraphIT	0.234	0.990 0.749	4.218
	(16,4)	GraphEM	0.325	0.928 0.302	4.185
		MLEM	0.792	0.016 0.031	4.741
		GraphIT	0.257	0.987 0.808	4.005
	(16,8)	GraphEM	0.322	0.907 0.403	3.916
		MLEM	0.650	0.031 0.061	3.964
		GraphIT	0.350	0.959 0.606	4.393
	(16,16)	GraphEM	0.362	0.905 0.494	5.572
		MLEM	0.710	0.063 0.118	3.129


est destinée au dépôt

Process of obtention of titanium nitride by rapid thermal annealing

Victor Process of obtention of titanium nitride by rapid thermal annealing. The titanium nitride (TiN) is an interesting material for plasmonic dispositive and metasurfaces. It can replace gold and silver in visible and near infrared spectral regions. Furthermore, due to its high thermal, chemical and mechanical stability the titanium nitride is interesting for applications where gold and silver are not stable.

The TiN can be obtained thanks to a lot of different physical techniques such as the physical vapor deposition (PVD) or the ammonolysis of titanium oxide (TiO 2 ). Regarding the ammonolysis, the TiO 2 layer need to react with ammonia gas at high temperature (more than 800°C) for several hours which hinder its industrial deployment.

In this presentation, the authors will present a promising alternative for the elaboration of titanium nitride thin films using a rapid thermal nitridation process (RTN) (figure 1). This technique consists in the substitution of oxygen of the TiO 2 by the addition of nitrogen under infrared flash of halogen light in a chamber containing ammonia [START_REF] Valour | Micro-Nanostructured TiN Thin Film: Synthesis from a Photo-Patternable TiO2 Sol-Gel Coating and Rapid Thermal Nitridation[END_REF][START_REF] Valour | Electrical and Mechanical Properties of TiN Thin Film Obtained from a TiO2 Sol-Gel Coating and Rapid Thermal Nitridation[END_REF]. This process is quick and allows to consider the industrial deployment. Nitridated layers have been characterized by X-ray diffraction, Raman spectroscopy, Transmission electronic microscopy (TEM and STEM), Electron Energy Loss spectroscopy (EELS). The optical and electrical properties will be presented as well. A plasmonic application will be presented.

  Vallejo-Otero, Arnaud Valour, Stéphanie Reynaud, Christophe Donnet, Émilie Gamet, Nicolas Crespo-Monteiro, Yves Jourlin To cite this version: Victor Vallejo-Otero, Arnaud Valour, Stéphanie Reynaud, Christophe Donnet, Émilie Gamet, et al.. Process of obtention of titanium nitride by rapid thermal annealing. ISNT : 11th International Symposium on Nitrides, May 2023, Saint-Malo (France), France. hal-04094503

Figure 1 .

 1 Figure 1. Illustration of a process of nitridation by RTN

  V. Vallejo-Otero1 , A. Valour 1 , S. Reynaud 1 , C. Donnet 1 , E. Gamet 1 , N. Crespo-Monteiro 1 et Y. Jourlin 1

1

Hubert Curien Laboratory, University of Lyon, Saint-Etienne, FRANCE


Anis M Haddouche
email: mohamed-anis.haddouche@imt-atlantique.fr

Fahima Djelil
email: fahima.djelil@imt-atlantique.fr

Christian Hoffmann
email: christian.hoffmann@univ-grenoble-alpes.fr

Nadine Mandran
email: nadine.mandran@univ-grenoble-alpes.fr

Cédric D'ham
email: cedric.dham@univ-grenoble-alpes.fr

Proposal of Indicators for Measuring Collaborative Writing in a Digital Learning Environment

Keywords: Indicators, Collaborative Writing, Computer-Supported Collaborative Learning (CSCL), Learning Analytics

Collaborative Writing (CW) is a common activity in education, which is being enhanced by the use of digital learning environments, leading to a growing research field in Computer-Supported Collaborative Learning (CSCL). In order to help teachers to monitor students CW, we propose two indicators that provide measures of student contributions to a text writing, namely balance of contribution and co-writing. We also identified CW strategies that are well defined in the literature. Moreover, we conducted a questionnaire evaluation to verify the interpretation of the indicators and the strategies by teachers in higher education context, using student reports edited in a collaborative digital environment called LabNbook, during physics and chemistry courses in undergraduate level. Results showed that teachers have a good interpretation of the indicators and strategies. This work contributes to research insights in CW, and motivates future work to design meaningful learning indicators.

INTRODUCTION

Providing students with group work is a common activity in education. Assigning group work is largely emphasized as a method that enables students to develop collaborative skills [START_REF] Sun | A tale of two communication tools: Discussion-forum and mobile instant-messaging apps in collaborative learning[END_REF]. Moreover, this may be enhanced by the use of digital technologies, providing students with several tools such as collaborative writing platforms [START_REF] Zhang | Assessing collaborative writing in the digital age: An exploratory study[END_REF]. This led to growing research in the field of Computer-Supported Collaborative Learning (CSCL) [START_REF] Chen | The role of collaboration, computer use, learning environments, and supporting strategies in cscl: A meta-analysis[END_REF]. More particularly, Collaborative Writing (CW) has gained an increasing research interest in the last few years [START_REF] Zhang | Computermediated collaborative writing in l2 classrooms: A systematic review[END_REF]. Students' interactions can be efficiently captured and stored, leading to a more fine-grained analysis of the CW process.

A large part of existing research aims at capturing student collaboration dynamics using mixed ap-dergraduate level. We address the following research questions: RQ 1.) What are the indicator metrics that allow to measure students' CW? RQ 2.) How can we deduce CW strategies from these indicators? RQ 3.) To what extent these indicators and strategies are interpretable by teachers? (How close is the relationship between these measures and their interpretation).

This paper is organized as follows. Section 2 provides a state of the art of existing research on CW. Section 3 introduces the proposed indicators and the co-writing strategies. Section 4 describes our research method. Section 5 details the calculations of the proposed indicators and their relationship to the cowriting strategies. Results are discussed in Section 6, conclusion and implications are derived in Section 7.

COLLABORATIVE WRITING

Collaborative Writing (CW) can refer to the production of a text by two or more writers (co-authoring) [START_REF] Storch | Collaborative Writing in L2 Classrooms. Multilingual Matters[END_REF]. More specifically, it is defined as a process involving substantive interactions between learners sharing decision making and responsibilities for a single produced document [START_REF] Zhang | Assessing collaborative writing in the digital age: An exploratory study[END_REF].

There is a long-standing interest in CW and different authors proposed early models and taxonomies of CW [START_REF] Posner | How people write together (groupware)[END_REF][START_REF] Lowry | Building a taxonomy and nomenclature of collaborative writ-ing to improve interdisciplinary research and practice[END_REF][START_REF] Storch | Collaborative Writing in L2 Classrooms. Multilingual Matters[END_REF]. With the fast-growing use of Online Learning Environments (OLE) and the ability to analyze data logs, several frameworks have been applied to examine behaviors, patterns and strategies of CW in different educational domains [START_REF] Onrubia | Strategies for collaborative writing and phases of knowledge construction in cscl environments[END_REF][START_REF] Sundgren | Visualizing online collaborative writing strategies in higher education group assignments[END_REF][START_REF] Olson | How people write together now: Beginning the investigation with advanced undergraduates in a project course[END_REF]. For instance, in the domain of second language (L2), a recent systematic literature review [START_REF] Zhang | Computermediated collaborative writing in l2 classrooms: A systematic review[END_REF] has examined more than one hundred studies revealing a strong research interest in the field of CW during the past decade.

One popular model in the literature that was proposed, is a dyadic interaction model with two constructs, namely equality, reflecting the learner's level of contribution and control over the task, and mutuality, reflecting the learner's level of engagement with each other's contribution [START_REF] Storch | Collaborative Writing in L2 Classrooms. Multilingual Matters[END_REF].

These two constructs also appear in less recent literature. For example, [START_REF] Dillenbourg | What do you mean by collaborative learning?[END_REF] compares concepts from the field of Human-Computer Collaborative Learning Systems (HCCLS), where an artificial agent collaborates with the human learner, and CSCL systems, where the computer supports collaboration between two human users. He argues that collaboration implies negotiation and emphasizes on the degree of symmetry in interactions between peers. He defines Symmetry of action as the extent to which the same range of actions is allowed to each agent [START_REF] Dillenbourg | Negotiation spaces in human-computer collaborative learning[END_REF] (the agent may refer to an artificial agent or a human peer). The term symmetry is borrowed from HCCLS domain to qualify an equilibrated balance of control implied by collaboration between the system and its user [START_REF] Dillenbourg | Negotiation spaces in human-computer collaborative learning[END_REF].

Mutual refinement is a second main characteristic defining collaboration through negotiation [START_REF] Baker | A model for negotiation in teachinglearning dialogues[END_REF]. This refers to specific strategies that exist for achieving agreement in the interaction (each agent successively refines the contribution of the other) [START_REF] Dillenbourg | Negotiation spaces in human-computer collaborative learning[END_REF]. This can be reflected in ways of editing text written by others and ways of coping with others editing one's own text [START_REF] Larsen-Ledet | Territorial functioning in collaborative writing[END_REF], or the degree of engaging with each other's ideas and each other's texts and providing scaffolding in producing joint writing [START_REF] Li | Explaining dynamic interactions in wiki-based collaborative writing[END_REF].

This led us to say that two key indicators may help to characterize collaborative work. The first one evaluates the symmetry in contributions, while the second one measures the degree of interactions between students.

INDICATORS AND STRATEGIES

We define balance of contribution as a metric that indicates how the student contributions are equal, well balanced or imbalanced. It is aligned on the one proposed by [START_REF] Olson | How people write together now: Beginning the investigation with advanced undergraduates in a project course[END_REF], called balance of participation which is a team contribution measure that reflects whether individuals' contributions are equal or imbalanced. As a balance metric, they considered one minus the variance of team members proportions in the collaborative work. The authors argued on variance as preferable over other ways of calculation such as Gini coefficient or Blau's index, for its simplicity and ease of interpretation. We based the calculation of the balance of contribution indicator on a variance metric, that provides a distribution between the authors' average contributions in terms of number of words in a final document. Balance of contribution indicator measures a student contribution at a word level, and co-writing indicator measures a student contribution at a sentence level. Indeed, balance of contribution measures only division of labor (number of words a student contribute to the text), and co-writing goes beyond that by measuring student co-construction of the text (modifying and adding words in sentences written by others in the text). The values of both indicators range between 0 and 1 (see section 5).

A possible exploitation of the indicators is the identification of CW strategies. We consider two strategies among the five proposed by [START_REF] Onrubia | Strategies for collaborative writing and phases of knowledge construction in cscl environments[END_REF], namely: 1) sequential summative text construction, i.e., one group member presents a document that constitutes an initial, partial or complete, proposal for the task resolution and the rest of the participants successively add their contributions to this initial document, without modifying what has been previously written, hence, systematically accepting what is added by other co-authors; 2) sequential integrative text construction, i.e., one group member presents a document that constitutes an initial, partial or complete task proposal, and the other group members successively contribute to this initial document, proposing justified modifications or discussing whether they agree with what has been previously written or not.

We use the term summative to appoint to a strategy where each student adds his text without modifying the text of the others where the result is being a juxtaposition of the individual contributions and, the term integrative, to appoint to a strategy where one student proposes an initial version and the other students contribute successively making modifications to the existing text [START_REF] Hoffmann | Development of actionable insights for regulating students' collaborative writing of scientific texts[END_REF].

The first is characterized by an explicit division of work between the team members, and the second by a co-construction of the text. Students may not necessarily follow one strategy but a mix of them.

RESEARCH METHOD

Our objective is, on the one hand, to construct metrics of students' CW process, and on the other hand, to verify that these metrics are intelligible for teachers. To conduct our research, we choose the design based research [START_REF] Wang | Design-based research and technology-enhanced learning environments[END_REF] and the associated guides [START_REF] Mandran | Thedre's framework: Empowering phd candidates to efficiently implement design-based research[END_REF]. This method proposes to build knowledge and associated tools in an iterative way by integrating all the actors in the field. Each iteration consists of developing the tool and the knowledge on the basis of the reality in the field and the opinion of the stakeholders.

Context

We use text documents produced by students on the LabNbook platform1 , which is a digital environment for learning experimental sciences in secondary and higher education (d'Ham et al., 2019). It provides useful tools for writing collaborative scientific documents (text and equations, drawing, data processing, etc.) and allows students to interact with each other (team building, shared workspace, internal messaging and chat, etc.). LabNbook operates in a locked co-editing mode [START_REF] Wang | Why users do not want to write together when they are writing together: Users' rationales for today's collaborative writing practices[END_REF], i.e., students can work simultaneously in the shared workspace but each document, called LabDoc, composing it can be edited only by one student at a time. For instance, a collaborative report may be composed of different sections, and each section corresponds to a LabDoc.

Data Collection

We use a questionnaire to measure to what extent the strategies and indicators are interpretable by teachers (RQ3). Our objective is to verify whether the definitions make sense to them.

We selected 12 LabDocs edited on the LabNbook platform by groups of students during physics and chemistry courses (undergraduate level) in Grenoble Alpes University (France). To distinguish student contributions, we highlight in each LabDoc with a color-code the text written by each student. We choose LabDocs uniformly with different writing strategies (integrative, summative, and mixed strategies). Respondent teachers are a total of 15, situated in IMT Atlantique in Brest (France), from Computer Science discipline and in Grenoble Alpes University in Grenoble (France), from Experimental Sciences disciplines. After providing the teachers with verbatim definitions of the strategies and of the indicators, we ask them to indicate for each LabDoc:

(1) the strategy they perceive (Entirely Summative (ES), Rather Summative (RS), Between Summative and Integrative (BSI), Rather Integrative (RI), Entirely Integrative (EI) and I don't know (DN));

(2) an estimate of the level of each indicator (Low (L), Medium (M) and High (H));

(3) a numerical value for each indicator between 0 and 1.

INDICATORS CONSTRUCTION

Text Sequences Matching Method

In order to study the evolution of a text document and compare pairs of text sequences, we use a Sequence Matcher method which has its origins in an algorithm published in the late 1980's by Ratcliff and Metzener under the name Gestalt Pattern Matching [START_REF] Ratcliff | Patternmatching-the gestalt approach[END_REF] 2 .

To illustrate this method, consider a text cowritten by two students A and B. Student A writes the first version which is then modified by student B. In order to qualify the evolution of the initial text to its final version, the approach consists firstly in finding the longest, in terms of number of characters, contiguous matching sequence (a set of words) that contains no useless elements, such as blank lines or white space. The same operation is then applied recursively for the sequences to the left and to the right of this longest contiguous matching sequence. Then, in order to qualify changes in the text, each sequence is tagged as follows: Equal (sequences are equal), Insert (sequence is inserted), Delete (sequence is deleted), Replace (sequence is replaced). For instance, let consider the following pair of text sequences written sequentially by two students:

• Student A: "LabNbook is a digital platform used by over 3500 students"

• Student B: "LabNbook is a platform used by more than 3500 students in France"

Student B contributes to the text after student A by inserting, deleting and modifying text. Table 1 gives the correspondent tag for each altered sequence.

Contribution Matrix

To quantify and track the evolution of the text, we use a matrix that gives for each word, the level of contribution of each of the students co-writing the text. We call this matrix, a contribution matrix, since it gives the student levels of contributions to a text. In this matrix, rows represent contributing students and columns represent words constituting the final text.

In order to provide a formal definition of this contribution matrix, let consider a text composed of N sentences, M words and co-written by K authors. Let also x i, j,l ∈ [0, 1], be the contribution level of an author i to a word l of a sentence j where i ∈ [1, K], j ∈ [1, N] and l ∈ [1, n j ]. Here, n j denotes the number of words in a sentence j.

In this K × M matrix, the level of contribution of an author is expressed by a score ranged between 0 and 1. It is set to 1 when the author writes entirely a word, and to 0, when the author doesn't contribute to a word. The total authors' contributions to a word is equal to 1.

Since our scoring method is based on the sequence matching method that identifies the longest contiguous matching sequence, deleted words do not result in a score, and replaced words result in a score that returns the ratio of contributions in terms of number of characters.

In the previous example, 2 students contribute to a text which is composed of 12 words ("LabNbook is a platform used by more than 3500 students in France"). Consequently, the contribution matrix gives scores ranged in 2 rows and 12 columns, as follows 1 1 1 1 1 1 0.3 0.3 1 1 0 0 0 0 0 0 0 0 0.7 0.7 0 0 1 1 , where the first row gives the scores corresponding to the contribution levels of the student A and the second row those of the student B.

Balance of Contribution Indicator

The balance of contribution indicator measures student's contribution at a word level. More precisely, it measures the degree to which students' contributions are equal (well-balanced) or not. The closer the value of this indicator is to 1, the more students' contributions are equal. At the opposite, the more it is closer to 0, the less the students' contributions are equal. This indicator is based on the variance in (3) of the students' average contributions, which reflects the distance between each student average contribution and the mean of all students' average contributions in (2). Note that, if the students contribute in a well-balanced way, then their scores are close to the mean and they contribute quite equally. Moreover, in order to penalize the case where the whole text is written by one student, this variance is normalized by (4). We define the balance of contribution indicator as follows

e(X) = 1 - K K -1 K ∑ i=1 x i, , - 1 K 2 (1)
where, for i = 1, . . . , K, x i, , is the average contribution of the student i to all words of the text. Note that, ∑ K i=1 x i, , = 1 and the mean of all students' average contributions is

x , , = 1 K K ∑ i=1 x i, , = 1 K . (2)
Indeed, as we consider that students contribute in a balanced way if their average contributions equals 1/K, we use in the construction of the balance of contribution indicator the variance

1 K K ∑ i=1 x i, , - 1 K 2 .
(3) However, in the case a student i writes alone the whole text, only one x i, , equals 1. Then (3) becomes

K -1 K 2 . ( 4
)
Finally, as mentioned above, in order to penalize the case where one student writes the whole text, we normalize the variance (3) with (4), which gives

K K -1 K ∑ i=1 x i, , - 1 K 2 .
This dispersion reaches its minimum value 0 when all students contribute equally (or in a balanced way) to the text, that is, where the average contribution of each student x i, , equals 1/K . It reaches its maximum value 1 when a student writes alone the whole text. In this case, only the average contribution of the student i, that is x i, , , equals 1. Therefore, in order to make this measure easier to interpret for teachers, we compute one minus this dispersion, which gives the balance of contribution indicator in (1).

Co-writing Indicator

Similarly to the balance of contribution indicator, which measures students' contributions at a word level, the co-writing indicator measures students' contributions at a sentence level (6). Thus, this indicator requires splitting the text into sentences. To this end, we use a rule-based (heuristic) approach for sentence segmentation, for its implementation simplicity [START_REF] Sadvilkar | PySBD: Pragmatic sentence boundary disambiguation[END_REF]. We use a sentence boundary detection tool based on a set of rules called Golden Rule Set 3 , hand-designed to determine sentence boundaries, such as punctuation.

3 https://github.com/diasks2/pragmatic segmenter When all sentences are written by a single student, this indicator equals 0 and it equals 1 when all sentences are co-written, in a balanced way, by all students.

The co-writing indicator is based on the balance of contribution indicator, but on a sentence level. According to the contribution matrix, we compute for each student his average contribution in all sentences. Thus, let x i, j, be, the average contribution of student i to the sentence j. Then, we obtain a K × N matrix of average contributions on a sentence level. For each sentence, we compute the variance of all authors' average contributions.

1 K K ∑ i=1 x i, j, - 1 K 2 . ( 5
)
Similarly to balance of contribution, we penalize the case where only one student writes a sentence alone, by computing the ratio between ( 5) and ( 4), and we calculate one minus the resulting quantity which gives the following dispersion measure

e j (X) = 1 - K K -1 K ∑ i=1 x i, j, - 1 K 2 .
Therefore, the co-writing indicator is given by

c(X) = N ∑ j=1 p j e j (X) where p j = n j M (6)
is the weight of the sentence j.

Indicators Property

The indicators present an inequality property. The cowriting indicator c(X) in ( 6) is lower or equal than the balance of contribution indicator e(X) in (1), i.e., c(X)e(X)

= K K -1 K ∑ i=1 x i, , - 1 K 2 - (7) N ∑ j=1 p j x i, j, - 1 K 2 ≤ 0 .
Indeed, equation ( 7) is non-positive as soon as

(x i, , ) 2 - n j ∑ l=1 p j (x i, j, ) 2 ≤ 0 , (8)
we have

(x i, , ) 2 = N ∑ j=1 p j x i, j, 2 ≤ N ∑ j=1 (p j x i, j, ) 2 . ( 9
)
Thanks to (9), an upper bound for (8) is given by

(x i, , ) 2 - n j ∑ l=1 p j (x i, j, ) 2 ≤ n j ∑ l=1 p j (p j -1) (x i, j, ) 2 ,
which is non-positive since p j ≤ 1 .

Collaborative Writing Strategies

As values of the indicators balance of contribution and co-writing are continuous and ranged between 0 and 1, we can represent co-written text documents in a two-dimensional (2D) plane plot (x-axis for balance of contribution and y-axis for co-writing), see Figure 1. Moreover, as co-writing indicator is always lower or equal than balance of contribution ((5.5)), the 2D plane is reduced to a triangle representing the area of potential plots for documents.

It is then possible to distinguish documents which are written in summative strategies (ES, RS) or integrative strategies (EI, RI) from documents written in a mixed strategy (BSI), by considering their location in the 2D plane. When the co-writing is low the strategy is summative. Therefore, documents written in summative strategies are plotted near the x-axis. When the strategies are integrative, the documents are plotted near the 2D plane diagonal and the balance of contribution and co-writing values are close. Documents written in a mixed strategy (between summative and integrative) are plotted in between.

RESULTS AND DISCUSSION

After providing the mathematical demonstration of the indicators (RQ1), and the illustration of the relationship between the indicators and the strategies (RQ2), the questionnaire results allow to evaluate the teacher interpretations of the strategies and the indicators (RQ3).

We compare the indicator values and levels estimated by teachers to the values computed using formulas (1) and ( 6 Teachers first indicate the perceived strategies. Ta-ble 2 compares the number of teacher responses received for each LabDoc and for each strategy. We report that 9/12 LabDocs are well classified with a majority of good responses (8 to 14). LabDocs 1 and 6 received respectively 6 and 5 well classifications while the majority of responses are spread-out on the other strategies. This is probably due to the difficulty for teachers to perceive a mixed strategy such as in LabDoc 6. LabDoc 1 comprises mathematical formulas, and this may influence the teacher perception.

For LabDoc 12, teachers answers are more spread out comparatively to LabDocs 1 and 6, but with 5 good answers and 1 answer as "I don't know". This can be explained by the fact that this LabDoc presents a very low collaboration (mainly written by a single student, a second contributor modified some words at its end) leading to a misinterpretation of the strategy by teachers. Secondly, teachers estimate for each LabDoc the level and value of the two indicators. Table 3 provides for each LabDoc the number of students having contributed to its writing, the machine value and machine level of each indicator, and the number of teacher responses for each indicator level.

We also calculate the distance between the machine value (MV) of an indicator θ ∈ [0, 1] and teachers' estimate, using the Root Mean Standard Deviation (RMSD) given by

RMSD = ∑ n 1 ( θi -θ) 2 n
where, for i = 1, . . . , n, θi ∈ [0, 1] is a teacher's estimate and n is the number of teachers who respond to the questionnaire. Concerning indicator levels, (8/12) of LabDocs are well classified by teachers regarding balance of contribution and (10/12) regarding co-writing. We also report that, for both indicators, mostly all miss-classified LabDocs have a value of RMSD greater than 0.3. This shows that when teachers succeed to perceive the correct indicator level for a LabDoc, they estimate well the indicator value. Moreover, regarding balance of contribution, most LabDocs for which the indicator levels are not correctly perceived, are those written by more than two students. Regarding co-writing, levels are not correctly perceived in LabDoc 1, this is due to mathematical formulas in this LabDoc.

From these results, and as an answer for RQ3, we can say that in some extent, the indicators and the strategies are interpretable by teachers and succeed to measure CW. Teacher misinterpretations are due to the specificity of some LabDocs, comprising mathematical formulas, or written in mixed strategies or of a very low collaboration. In fact, teachers examined manually the LabDocs, and these specific characteristics led to misconceptions.

CONCLUSION AND IMPLICATIONS

In this work we contribute to research in CW with new indicators that provide measures of student contributions in a text writing. A first contribution consists in two indicators that allow to measure how student contributions are equal (balance of contribution), and how students interact with each other text contributions (co-writing).

A second contribution consists of identifying CW strategies, namely summative and integrative strategies that can be useful in pedagogical settings, for example to distinguish between cooperative (summa-tive) and collaborative (integrative) works [START_REF] Onrubia | Strategies for collaborative writing and phases of knowledge construction in cscl environments[END_REF]. We showed that these writing strategies can be derived from the two proposed indicators.

A third contribution is the evaluation of the interpretation of the indicators and strategies, allowing to compare teacher perceptions of the indicators with machine calculations, and to verify the perceived strategies. Results showed that teachers have a good interpretation of the indicators and the strategies. In fact, teachers were mostly able to estimate correctly the indicators based on their textual definitions, as well as the writing strategies. We conclude that the indicators may help for monitoring student CW, and characterizing the degree and strategies of collaborative work.

This work is of scholarly and practical implications. It provides interesting insights on the design of learning analytics tools that allow for a meaningful reporting of CW dynamics between peers. This may strengthen formative evaluation and provide learners with quick feedback during their collaborative work and reports co-writing. This work is not without limitations. Our approach can produce misleading results when collaborative documents are composed of mathematical formulas. More advanced text mining methods can improve formulas detection. Moreover, it would be relevant to collect more data with teachers and students to demonstrate the value of the proposed indicators, and investigate better insights into CW and collaborative work.

Future work is motivated to measure student cognitive contributions, by considering the meaning of sentences composing the text. Our approach can be improved by considering semantic sentence detection using machine learning models. This may better inform about the collaboration process and go beyond contribution measures in terms of number of words added or modified in a text. Ultimately, it would be worthwhile to study how to enable teachers to take actions from useful visualizations of indicators and strategies, and evaluate their acceptability, usability and utility.

Figure 1 :

 1 Figure 1: A 2D plane illustrating the relationship between the indicators and the CW strategies.

  ). Indicator levels are set among Low for indicator values in ]0, 1/3], Medium for values in [1/3, 2/3[ and High for values in [2/3, 1].

Table 1 :

 1 Illustration of the tagging operation.

	Tag	Student A	Student B
	Equal	"LabNbook is a"	"LabNbook is a"
	Delete	"digital"	" "
	Equal	"platform used by" "platform used by"
	Replace "over"	"more than"
	Equal	"3500 students"	"3500 students"
	Insert	" "	"in France"

Table 2 :

 2 Teacher responses regarding the strategies they perceived. In bold, LabDocs that are misclassified by the majority of teachers. The wright strategy is ST.

	LabDoc ST ES RS BSI RI EI DK
	1	EI	0	2	4	3	6	0
	2	ES	11	2	1	1	0	0
	3	EI	0	1	1	1 12	0
	4	BSI	0	1	8	5	1	0
	5	EI	2	0	0	2 11	0
	6	BSI	0	8	5	2	0	0
	7	ES	12	1	0	0	2	0
	8	RI	0	0	4	7	4	0
	9	RS	0	8	5	1	0	1
	10	ES	14	0	1	0	0	0
	11	ES	14	1	0	0	0	0
	12	EI	1	4	2	2	5	1

Table 3 :

 3 Comparing teacher perceptions and estimations of indicator levels and values with the machine values (MV) and levels (ML). In bold, LabDocs for which a majority of answers are wrong. "NB" is the Number of Students (co-writers).

				Balance of contribution			Co-writing
	LabDoc NB MV RMSD ML L M H MV RMSD ML L M H
	1	2	0.49	0.32	M	0	7	8 0.41	0.42	M	0	4 11
	2	2	0.87	0.19	H	0	6	9 0.11	0.15	L	13 1	1
	3	2	0.84	0.24	H	1	2 12 0.76	0.28	H	1	0 14
	4	2	0.99	0.21	H	0	1 14 0.41	0.22	M	3	9	3
	5	3	0.49	0.36	M	15 0	0 0.44	0.35	M	12 1	2
	6	2	0.99	0.19	H	0	1 14 0.50	0.2	M	7	7	1
	7	2	0.51	0.3	M	13 2	0 0.07	0.09	L	14 1	0
	8	2	0.99	0.3	H	0	7	8 0.74	0.21	H	0	6	9
	9	4	0.48	0.27	M	13 2	0 0.17	0.12	L	13 2	0
	10	3	0.92	0.39	H	1 11 3 0.09	0.1	L	13 2	0
	11	2	0.29	0.19	L	15 0	0 0.04	0.04	L	15 0	0
	12	2	0.14	0.07	L	15 0	0 0.09	0.09	L	14 1	0
	RMSD take values in [0.07, 0.39] regarding balance of
	contribution and in [0.04, 0.42] regarding co-writing.

https://LabNbook.fr/

For implementation, we use the Python library Difflib https://github.com/python/cpython/blob/main/Lib/difflib.py

ACKNOWLEDGMENTS

The project is co-financed by the Brittany Region within the research program SAD and the Finistère Department Council within the research program APRE. We thank all the teachers having participated to the study.