Milan Courcoux-Caro

Charles Vanwynsberghe

Cédric Herzet

Alexandre Baussard

Sequential sensor selection for the localization of acoustic sources by Sparse Bayesian Learning

This paper deals with the design of sensor arrays in the context involving the localization of a few number of acoustic sources. Sparse approximation is known to be eective to nd the source locations, but it depends on dierent array characteristics such as the number of sensors and the array geometry. The present paper tackles this array design problem under the form of a sequential sensor selection procedure. The proposed method alternates between two steps. One step involves a source localization estimator, given a current set of measurement points, to obtain the estimation variance. Then the other step selects the new point where a future measurement will maximally decrease the variance from previous step.

As such, the procedure can be applied online. Both numerical and experimental studies are conducted in an indoor neareld conguration. Results show that the proposed approach performs better than oine state-of-the art methods, and the presented empirical study reveals a better robustness to the model mismatches originating from the room reections.

I. INTRODUCTION

Source localization from a limited number of pointwise measurements is a common inverse problem in acoustic or radar applications. The literature exhibits a plethora of techniques to solve this problem.

To name a few, it goes from methods of beamforming 9,41 , deconvolution [START_REF] Brooks | A deconvolution approach for the mapping of acoustic sources (DAMAS) determined from phased microphone arrays[END_REF] , to estimators relying on the help of a prior through a Bayesian framework 2 , or Compressed Sensing (CS) when there are few sources in the eld [START_REF] Xenaki | Compressive beamforming[END_REF] .

However, irrespective of the chosen technique, the quality of source localization is impacted by the conguration of the array. To improve the reconstruction accuracy, the classical antenna theory suggests the trivial solution of increasing the number of sensors, placed at randomized positions to avoid a regular geometry structure. Such an approach remains valid, especially when the chosen method is sensitive to side and grating lobes such as beamforming. However, in this context the randomization of sensor positions does not meet any form of optimality. Moreover, depending on the context, the number (resp. positions) of sensors may be constrained by cost (resp. structural) limitations.

a) This work was also done while the author was aliated at (a) An extensive literature depicts how to optimize the array geometry. One well-known approach relies on a renement of the directivity pattern, such as the reduction of the side lobe levels [START_REF] Boeringer | Particle swarm optimization versus genetic algorithms for phased array synthesis[END_REF][START_REF] Yan | Sidelobe reduction in arraypattern synthesis using genetic algorithm[END_REF] . The geometry is obtained as a solution of a non-convex minimization problem, typically via a genetic algorithm [START_REF] Haupt | Thinned arrays using genetic algorithms[END_REF][START_REF] Courtois | Optimisation par algorithme génétique de la géométrie d'antenne pour la localisation de sources[END_REF][START_REF] Liu | Optimal sensor placement for spatial lattice structure based on genetic algorithms[END_REF] , swarm optimization [START_REF] Khodier | Linear array geometry synthesis with minimum sidelobe level and null control using particle swarm optimization[END_REF] or simulated annealing [START_REF] Murino | Synthesis of unequally spaced arrays by simulated annealing[END_REF][START_REF] Trucco | Thinning and weighting of large planar arrays by simulated annealing[END_REF] . This approach is popular and remains relevant to improve beamforming eciency, since the latter is inherently related to the beampattern shape.

However, it is not relevant for other methods than beamforming. For example in CS, the quality of sparse source localization depends on the incoherence of the acoustic propagation operator. In this case, the geometry is preferably designed based on a statistical approximation of the restricted isometry property, or the null space property [START_REF] Gaumond | Application of statistical reduced isometry property to design of line arrays for compressive beamforming[END_REF] .

Finding an optimal array geometry can be tackled by ideas stemming from experimental design 13,[START_REF] Pukelsheim | Optimal design of experiments[END_REF] , and is known under the name of sensor selection in the literature. The goal is to minimize the variance of some estimator for a given number of sensors and a collection of candidate positions. Although attractive, acoustic array design by sensor selection seems less common than beampattern-based methods. One example, in [START_REF] Chepuri | Sparsity-promoting sensor selection for non-linear measurement models[END_REF] , illustrates an optimal array deployment for localization by time of arrival measurements. In 18 , the number of measurements is sequentially augmented by moving a microphone array to a new optimal position, and the source reconstruction is rened at each iteration. Despite being elegant, this approach is specic to spread sources due to the use of Gaussian prior.

The aim of the present work is to introduce a sequential sensor selection when there are few acoustic sources.

It is derived from the Bayesian framework and relies on two key points. First it computes a solution of the source vector by a sparsity-enhancing hierarchical model [START_REF] Yang | O-grid direction of arrival estimation using sparse bayesian inference[END_REF] . Second, it selects the next sensor position which maximally reduces the variance of the previously computed solution.

With such a sequential procedure, the optimal selection is informed by the previous measurements [START_REF] Ji | Bayesian compressive sensing[END_REF] , and successful applications were demonstrated such as in MRI 36 or temperature monitoring [START_REF] Krause | Near-optimal sensor placements in gaussian processes: Theory, ecient algorithms and empirical studies[END_REF] . The present work investigates the eciency of the idea in acoustic source localization, with a validation in an experimental, neareld indoor scenario. To the best of authors' knowledge, this is the rst paper exhibiting in situ experiments of sensor selection applied to source localization with real data.

The paper is organized as follows. In sec. II, the theory of the proposed approach is described.

II. SENSOR SELECTION FOR SOURCE LOCALIZATION

A. Model description

Let the array of M sensors lie in the 3-dimensional space, with the m-th sensor position denoted p m ∈ R 3 . The set M = {p 1 , . . . , p m , . . . , p M } contains the sensor positions, meaning that |M| = M . The main goal of our work is to determine a subset K ⊆ M of K sensor positions leading to good sourcelocalization performance. We suppose that there are S sources at frequency f to localize, and we approximate the continuous source domain by a grid of N points.

Thus, we consider a linear on-grid model with additive noise to describe each measurement. For the m-th sensor the model writes:

y(p m ) = a(p m )s + w m , (1)
with s ∈ C N ×1 the vector of unknown source amplitudes to be estimated, a(p m ) ∈ C 1×N the row vector describing the propagation of the sources from the grid to the sensor position p m at frequency f , and w m ∈ N C(0, σ 2 ) is a random noise.

The acoustic eld is supposed to contain a few number of sources, with the condition that S ≪ N . Accordingly the signal s is S-sparse, i.e. it has S non-zero elements. Note that there is no particular constraint on the propagation model, excepted linearity: the row vector a(p m ) may depict conditions of far or near eld, with or without specic boundary conditions depending on the eld of application e.g. indoor with multipaths or underwater in shallow water [START_REF] Jensen | Modern Acoustics and Signal Processing Computational Ocean Acoustics[END_REF] .

Let K ⊆ M to be a subset of sensors from the array, meaning that |K| ≤ M . Dening the measurement vector as y K = y(p 1 ), . . . , y(p m ), . . . , y(p |K| ) T , the model (1) in matrix form reads:

y K = A K s + w K , (2)
with A K = a(p 1 ) T , . . . , a(p m ) T , . . . , a(p |K| ) T T the |K| × N forward propagation matrix, and w K = w 1 , . . . , w m , . . . , w |K| ) T the |K| × 1 noise vector. We assume that the elements of w K are independent and identically distributed:

w K ∼ N C(0 |K| , σ 2 I |K| ).
In this paper, we consider a Bayesian framework to devise our sensor-placement procedure.

To this end, we need to dene a prior model P(s) which properly accounts for the sparse nature of s. A typical choice for P(s) is the so-called Bernoulli-Gaussian model (see e.g. [START_REF] Herzet | Bayesian pursuit algorithms[END_REF][START_REF] Herzet | Sparse representation algorithms based on mean-eld approximations[END_REF][START_REF] Schniter | Fast Bayesian matching pursuit[END_REF][START_REF] Soussen | From BernoulliGaussian Deconvolution to Sparse Signal Restoration[END_REF][START_REF] Zayyani | Decoding real eld codes by an iterative Expectation-Maximization (EM) algorithm[END_REF]. One standard formulation of this model reads:

P(s) = N n=1 P(s n ), (3)
with

P(s n ) = bn∈{0,1} P(s n |b n )P(b n ). (4)
P(b n ) = Ber(p) is a Bernoulli distribution with param- eter p ∈ [0, 1] and P(s n |b n ) = CN (0, σ 2 bn
) is a zeromean complex Gaussian density with variance σ 2 bn . It can be easily seen that typical realizations of the model are (quasi-) sparse vectors when p ∼ 0, σ 2 0 ∼ 0 and σ 2 1 ≫ 0. Interestingly, several authors (see [START_REF] Herzet | Bayesian pursuit algorithms[END_REF][START_REF] Soussen | From BernoulliGaussian Deconvolution to Sparse Signal Restoration[END_REF] have emphasized the connection between a maximum a posteriori estimation problem involving a Bernoulli-Gaussian prior and the ℓ 0 -penalized problem, ubiquitous in the eld of sparse representations. The Bernoulli-Gaussian model can thus be seen as an ideal probabilistic modelling of the sparse nature of the source vector s.

In this work, we use the Bernoulli-Gaussian prior as a starting point to motivate our sensor-selection methodology. In the next sections, we will then consider some simplications to circumvent the computational bottleneck induced by this model.

B. The sensor selection problem

In the context of source localization, the sensors must be positioned in such a way that the uncertainty in the number, positions and amplitudes of the sources is minimized when the measurements are collected.

In the Bayesian framework considered in this paper, the uncertainty in these quantities is directly related to the spread of the posterior distribution p(s|y K ). For instance, if p(s|y K ) reduces to a Dirac's delta, s can be perfectly recovered from y K and there is no residual uncertainty on the source locations upon the observation of y K . Contrarily, when the mass of p(s|y K ) is spread over the whole domain C N , the uncertainty about s is high and the error made by any point-estimate sK is likely to be large. Since working with the full distribution p(s|y K ) is an intractable task, it is common to map the latter to some scalar-valued cost function, see 4,32 . In this paper, we use the following function:

f (K) = log det(Σ K ), (5)
where

Σ K ≜ E s|y K {(s -E s|y K {s})(s -E s|y K {s}) H } (6)
is the covariance matrix of s|y K . This gure of merit is used in the so-called D-optimal design procedures and is widely known in the literature [START_REF] Pukelsheim | Optimal design of experiments[END_REF] . There are two motivations for using this function. First the covariance matrix fully describes the estimator variance, so a selection based on this matrix will naturally control the estimator uncertainty. Second, the covariance matrix describes the shape of the condence ellipse in the N-dimensional space that contains the error estimation [START_REF] Joshi | Sensor selection via convex optimization[END_REF] . The goal of the function is to capture the ellipse shape with one scalar.

In D-optimality, the function f (K) provides the log volume of this ellipse.

Note that relying on the covariance only has a limitation: it does not take the estimator bias into account.

If the mean square error expression of the considered estimator exists, the optimal selection can be based on it so that both bias and variance are reduced e.g. see 33 . This direction is out of scope in this paper, because the mean square error is unknown. Thus, there is no guarantee of controlling the estimator bias via D-optimal design.

The main question addressed in the paper is as follows: which sensors in M should we select in order to minimize f (K) when |K| = K? This question formally write as follows:

Problem 1 Given a set of potential sensor positions M and a number of sensors K ≤ M , nd the subset of sen-sor positions K such that

K = arg min K⊆M f (K) subject to |K| = K. (7)
Unfortunately, solving this problem is an intractable task. Bottlenecks are located at two levels. First, ( 7) is a discrete-valued optimization problem and, in the current state of knowledge, there is no polynomial-time procedure able to solve it. Second, the evaluation of the cost function f (K) is intractable because the Bernoulli-Gaussian prior induces a summation over the 2 N possible congurations of the Bernoulli variable when p(s|y K ) has to be evaluated.

In the rest of this paper, we propose a methodology for sensor selection which circumvents these two bottlenecks by resorting to some approximations. In the next two subsections, we show that the posterior covariance matrix Σ K (and thus the cost function f (K)) can be efciently evaluated when the prior model P(s) is relaxed to some hierarchical Gamma-Gaussian model. In section III, we consider a greedy procedure to compute an approximate solution of (7).

Our approach can thus be seen as a tractable procedure to compute an approximation of the solution of (7).

C. Posterior covariance with Gaussian prior model

As a starting point, let us focus on the structure of the posterior covariance matrix when considering a

Gaussian prior model:

P(s|α) = CN (0, diag(α) -1 ), (8)
where α ∈ R N + denotes the inverse-variance of each ele- ment of s. ? In this case, it is well known that the posterior covariance Σ K takes the form 2 :

Σ K = σ -2 A H K A K + diag(α) -1 -1 . (9)
We note that Σ K is well-dened even when the problem is underdetermined (M < S) or ill-posed. Moreover, it has a nice closed-form solution which makes it attractive for our sensor-placement problem (7).

Unfortunately, it is well-known that Gaussian model (8) does not promote sparse realizations. It is therefore a poor prior to account for the sparse nature of the source vector s. In fact, the only way to enforce the n-th element of the source vector s to be close to zero with high probability is to set α -1 n ≃ 0. If α is user-provided, this means that the source locations should be known in advance to obtain a proper prior model on s.

In the next section, we investigate a solution to benet from the closed-form expression of the Gaussian posterior covariance matrix while accounting for the prior lack of knowledge of α.

D. Hierarchical Gamma-Gaussian prior model

In the previous section we emphasized that Gaussian prior models lead to desirable closed-form expressions of the posterior covariance matrix. Unfortunately, this prior is only of practical relevance when the parameter vector α is properly tuned to account for the sought source positions s. In this section, we circumvent this issue by: i) considering a Gaussian prior model on s as in (8); ii) including α as an unknown quantity in our estimation problem. The rationale for this strategy is as follows:

even if a proper choice for α is not known in advance, a relevant value for this parameter may be learned from the received observations.

More precisely, we consider the following steps. First, since we assume that α is an unknown quantity and we operate in a Bayesian context, we must dene a prior model on α. In this work, we rely on the well-known hierarchical model of Sparse Bayesian Learning to enforce sparsity [START_REF] Tipping | Sparse bayesian learning and the relevance vector machine[END_REF] , and dene the following prior:

P(α|ρ) = N n=1 Γ(α n |1, ρ), (10)
where Γ(α n |1, ρ) denotes a Gamma distribution with shape parameters 1 and ρ. We set ρ ≫ 1 which leads to a at distribution, making all the values of α a priori roughly equally probable see for example [START_REF] Yang | O-grid direction of arrival estimation using sparse bayesian inference[END_REF] . We note that combining ( 8) and ( 10) is tantamount to dening a Gaussian-Gamma hierarchical prior model on (s, α).

Here, the benet from such a prior is that it is noninformative.

Second, we compute a maximum a posteriori estimate of α, i.e, αK = arg max

α∈R + P(α|y K ). (11)
We note that problem (11) with a Gaussian-Gamma prior model has already been considered in the context of sparse representations (see [START_REF] Tipping | Sparse bayesian learning and the relevance vector machine[END_REF] ) and source localization [START_REF] Gerstoft | Multisnapshot sparse bayesian learning for DOA[END_REF][START_REF] Yang | O-grid direction of arrival estimation using sparse bayesian inference[END_REF] . In particular, it was emphasized that the solutions of ( 11) typically contain a few large components and many close-to-zero elements. This is in good agreement with our target problem where one searches to turn on only the elements of α corresponding to emitting sources. A numerical procedure, dubbed Sparse Bayesian Learning (SBL), instantiating the recursion of an expectation-maximization algorithm to problem (11) was also proposed in [START_REF] Tipping | Sparse bayesian learning and the relevance vector machine[END_REF] . We will consider this algorithm to compute αK and refer the reader to 45 for implementation details.

Third, we approximate the posterior covariance matrix by conditioning on α = αK :

Σ K ≃ E s|y K , αK {(s -E s|y K , αK {s})(s -E s|y K , αK {s}) H }.
Interestingly, if the value of α is given, the prior model becomes Gaussian and one recovers the expression of the posterior covariance matrix given in (9), that is

Σ K ≃ σ -2 A H K A K + diag( αK ) -1 -1 . ( 12
)
We note that the only dierence between ( 9) and (12) is that the unknown parameter α has been replaced by the posterior estimate αK . We therefore recover the nice analytical expression of the covariance matrix obtained in the case of a Gaussian prior model but with a value of the parameter α estimated from y K .

Towards a data-informed sensor selection

We remark that SBL is particularly attractive because the computed covariance is learned posterior to the sensor readings. Indeed, SBL helps to go towards a sensor selection strategy that is driven by both the model and the measured data, because αK is computed poste- rior to y K .

From this key point, the next section proposes a data-informed sensor selection method, based on an iterative procedure.

III. DATA-INFORMED SENSOR SELECTION: A GREEDY IT-ERATIVE ALGORITHM

Problem 1 is combinatorial and is therefore not tractable with a high-dimensional sensor set M. For applicability to large-scale arrays, we choose to follow a greedy heuristic approximating the solution. There are two motivations:

1. taken together the steps of selection have a low

complexity in O(KM N 2 ) justied below,
2. the framework is naturally online and datainformed, since each iteration allows to take the measurement from the new selected sensor into account.

Accordingly, the algorithm achieves a sequential design, and increments the size of K until |K| = K.

A. The one-sensor problem

Assume one wishes to add one extra sensor to a subset K ⊆ M. In the spirit of (7), it seems reasonable to select the sensor leading to the maximum decrease of cost function f , that is

p = arg max p∈M\K f (K) -f (K ∪ p). (13)
In this section, we show that using covariance-matrix surrogate (12) with an additional approximation leads to a tractable analytical solution to this problem.

Let the function f (K) be f (K) = log det(Σ K ) with ( 14)

Σ K = σ -2 A H K A K + diag( αK ) -1 -1 , ( 15
)
where αK is the solution of [START_REF] Demontis | Block-sparse approach for the identication of complex sound sources in a room[END_REF]. With these denitions, a brute-force evaluation of (13) turns out to be computationally intensive since it requires to evaluate αK∪p for any p ∈ M\K. To improve the computational eciency of our proposed method, we consider hereafter the following additional approximation:

∀p ∈ M\K : αK∪p ≈ αK . ( 16
)
The key point of this approximation is to access a tractable computation of the new covariance Σ K∪p with an hypothetical additional sensor, without the explicit need of the new measurement y(p) inside the expression.

As a result, this leads to

Σ -1 K∪p = σ 2 A H K∪p A K∪p + diag( αK∪p ) -1 ≈ σ 2 A H K∪p A K∪p + diag( αK ) -1 = σ 2 A H K A K + diag( αK ) -1 + σ-2 K a(p) H a(p) = Σ -1 K + σ 2 a(p) H a(p). (17)
Using the matrix determinant lemma [START_REF] Ding | Eigenvalues of rank-one updated matrices with some applications[END_REF] , we obtain

det(Σ K ) ≈ det(Σ K∪p )(1 + σ 2 a(p)Σ K a(p l ) H ). ( 18
)
Finally, plugging this expression into (13), we end up with the following simplied optimization problem:

p = arg max p∈M\K log 1 + σ 2 a(p)Σ K a(p) H . ( 19
)
Note that this sequential form of D-optimal design was derived by Wynn et al in [START_REF] Wynn | The sequential generation of d-optimum experimental designs[END_REF] in the case of linear models.

Through this form, computing the log volume reduction is ecient as it requires O(N 2 ) operations to score the function for each candidate (once Σ K has been evaluated).

B. Remarks and connections with related literature

Remark 1

The eq. ( 19) was also derived by information-theoretical tools [START_REF] Seeger | Compressed sensing and bayesian experimental design[END_REF] . It helps to understand its optimality as a choice of the most informative sensor. The dierential entropy, by denition, describes the entropy of a random variable admitting a distribution with respect to the Lebesgue measure. In the considered hierarchical model for SBL with measurements from K, the conditional probability density function P(s|y K , αK ) is a multivariate complex Gaussian distribution. In this case, the dierential entropy writes: Applied to the multivariate normal posterior distribution of s in SBL with measurements in K, it writes:

h K = - C N f (s) log f (s)ds, with f (s) = e (s-s K ) H Σ -1 K (s-s K ) π det(Σ K ) . ( 20
)
In this case, it is well known that 8 (Ch. 8):

h K ≡ log det(Σ K ), (21)
where ≡ means equality up to a constant. According to the same guideline, the dierential entropy with one extra sensor amounts to h K∪p ≡ log det(Σ K∪p ). As a result, the score in eq. ( 13) also amounts to the reduction of entropy h K -h K∪p , i.e. the information gained thanks to a sensor reading at p. With the help of approximation in eq. ( 16), the information gain nearly equals the eq. ( 19).

Remark 2

One geometric interpretation of the constrained log volume reduction is the following. Supposing that the row norm ∥a(p)∥ is invariant, we search the direction in C N colinear to the major radius of Σ K , under the constraint that a(p) belongs to the manifold structured by the physical model. Note that, according to [START_REF] Ji | Bayesian compressive sensing[END_REF] , if a lies in C 1×N without any physical constraint, (i) it is not function of p, and (ii) it maximizes the log volume reduction when it is equal to the eigenvector paired with the largest eigenvalue of Σ K .

C. Final algorithm

As a whole, the proposed method is an online greedy sensor selection. Based on the covariance a posteriori from SBL, it is naturally data-informed. Note that it inherently calls SBL, so the selection jointly performs an online source localization. Indeed, the source vector sK can be estimated based as a function of Σ K as:

sK = Σ K A H K y K . ( 22
)
The method is summarized in the algorithm 1.

Input : K = ∅, K, M, yM Note that an initialization step is necessary: K cannot be empty in the present case, since SBL needs data, including at the rst iteration. One initial sensor is enough in theory, but not enough to extract a spatial information on the source eld with SBL. Thus, two initial sensors are picked at random and included in K.

[a(p m )] n = 1 ||p m -r n || 2 e -ı 2πf c 0 ||pm-rn||2 , ( 23
)
where f is the frequency in Hz, c 0 the sound speed in m/s, and r n the position of the n-th source on the grid.

Obviously, such a model is mismatched indoor. It can be rened, by taking the presence of walls into account.

Reections from the 6 surfaces (i.e. the 4 walls, ceiling and ground) are integrated to the model with the method of image sources 1 . To do so, each point of the original source grid is duplicated with respect to the surface positions, as depicted in g. The image source method is convenient, because it is linear and its analytical expression is easily calculable.

Indeed, the general form of its related propagation model can be written as [START_REF] Allen | Image method for efciently simulating small-room acoustics[END_REF][START_REF] Ribeiro | Using reverberation to improve range and elevation discrimination for small array sound source localization[END_REF] :

[a(p m )] n = 1 ||p m -r n || 2 e -ı 2πf c 0 ||pm-rn||2 + J j=1 β j 1 ||p m -r j n || 2 e -ı 2πf c 0 ||pm-r j n ||2 , (24)
with J the total number of eective reections, r j n the position of the j-th image source corresponding to the original source point r n , and β j each eective attenuation coecient.

We assume that the absorption by the walls does not induce any phase shift, and does not depend on the wave incidence angle. So β j are real constant numbers for each index j. The reection coecients are roughly estimated by nding a high correlation between the measurements and synthesized data. We x β j = 0.75 for the ceiling and walls, and β j = 0.5 for the ground as coecients of the rst order. Other coecients of the z-th order are derived accordingly, by taking into account the appropriate reections occurring in the acoustic path. By awareness, we veried that small variations of β j lead in marginal impacts on the provided results in sec. VI.

The addressed model (24) will be used to approximate the indoor propagation but obviously the level of approximation will introduce model mismatches, the free-eld model ( 23) being the coarsest choice. Nevertheless, promising works show that SBL is robust to mismatches because of its adaptive nature that results from the underlying probabilistic model 15,16 . To follow up on this work, the experiments in sec. VI will analyze if the data-informed sensor selection leverages this adaptive nature to enhance the localization accuracy.

V. EVALUATION SETUP OF THE PERFORMANCE

In this section we briey introduce the state-of-theart methods used as references to evaluate the performances of the proposed approach. Before that, we dene the criterion used for the comparison.

A. Evaluation criterion

To compare the performances of the proposed method against the state-of-the art ones (described below), we propose to use the following criterion: once the antenna array has been selected for a given number of sensors, the SBL is used to provide the complex source vector sK [START_REF] Herzet | Sparse representation algorithms based on mean-eld approximations[END_REF] on the grid. Note that our approach in- herently involves the SBL computation online.

Starting from the prior knowledge that we have S sources, we select the S highest peaks from the absolute value of sK elements (the SBL doesn't necessarily gives spikes, that is why we look for the positions of the local maxima). Then, we compute the euclidean distance between the true source positions and these positions. If this distance is lower than some given radius threshold, then we consider that the estimated source location is a true positive. Let N tp be the number of true positives. Moreover, since we are going to achieve Monte Carlo simulations with N r realizations, the result will be displayed as empirical mean of true positives in percent:

True positive = 100 N r Nr i=1 N tp,i S [%], (25)
where i stands for the index of the Monte Carlo realization. In addition, the Jaccard index will be plotted because it also takes False positives and False negatives into account. We recall, from its expression

Random selection (SelRand)

Since our approach is based on the sparsity prior, it seems appropriate to consider solutions from the CS literature. In this framework, making random linear measurements of sparse vectors is known to be eective to obtain good reconstruction guarantees. However, in the present acoustic problem, the matrix A K is deterministic and its elements can not deliberately randomized. But studies show that the requirements from the CS framework can be approximated by selecting the sensor positions at random 6,43 . To compare this strategy with the proposed approach, we simply replace in Algorithm 1 the optimal selection p by a uniform random selection in the set M\K.

We still compute SBL at each iteration though, in order evaluate the approach according to the (increasing) number of selected sensors.

Beampattern optimization with genetic algorithm (SelGen)

Array design can be achieved by an optimization problem seeking to reduce most of the time the Maximum Sidelobe Level (MSL) of the beampatterns 20 .

Shortly, the goal is to maximize the dierence between the main lobe and the sidelobe levels. This kind of approach is especially relevant for source localization using beamforming. The optimization problem, which is non-convex, can be solved using several approaches. In this paper, we use the MSL-based beampattern optimization based on a genetic algorithm proposed in [START_REF] Courtois | Optimisation par algorithme génétique de la géométrie d'antenne pour la localisation de sources[END_REF] . Once the array is designed for a given number of sensors, we perform SBL and get the evaluation criterion (25).

Note that this method is oine, and not inherently iterative. So, to compare with the online sequential selections, the arrays are independently designed for dierent values of K. Finally, because the present experimental setup is dierent from the one in [START_REF] Courtois | Optimisation par algorithme génétique de la géométrie d'antenne pour la localisation de sources[END_REF] , some adjustments should be done regarding the initial choice of three parameters: the number of generations, of individuals, and the mutation coecient numerical values are given in section VI.

Convex relaxation with the Fisher information matrix (SelRelax)

The last method considered in this paper is the one proposed in [START_REF] Chepuri | Sparsity-promoting sensor selection for non-linear measurement models[END_REF] . The goal of this method is to select the positions of a given number of sensors by maximizing the smallest eigenvalue of the Fisher Information Matrix (FIM) via convex relaxation. Note that maximizing the smallest FIM eigenvalue is known as E-optimal design.

In This method is model-informed only and, more importantly, is not sequential de facto like the method in sec. V B 2. Thus, we proceed in the same way by designing the arrays independently for dierent values of K.

VI. RESULTS

In this section we present results exploiting both simulated and experimental data, in the scenario of source localization indoor in near eld, thoroughly described in sec.IV. Two main analyses are conducted: the eect of model mismatch and the performance evaluation of the proposed procedure against state-of-the-art approaches (c.f. section V B).

For the eect of model mismatch, we focus on how the steering matrix A K is dened regarding the reverberation. Based on the dened model in sec. IV B, we evaluate the dierence of performances of source localization indoor, when the model either is free eld (as eq. ( 23)), or takes reections into account (as eq. ( 24)).

Concerning the implementation of the oine algorithms SelGen and SelRelax, described in sections V B 2 and V B 3, arrays are designed for K starting from 50 with a step of 50, and stopping at 1000 and 650, respectively. Note that the latter is stopped early because the convex optimisation by CVX fails beyond. Finally, some parameters of the genetic algorithm (SelGen) from 28 are adjusted: the number of generations and individuals are higher (250 against 100) and the mutation coecient is lower (0.0001 against 0.01).

A. Results from numerical studies

We generate a set of N r = 300 realizations by randomly selecting S = 4 source positions of the 25 plotted in g.2(b). The plane containing the sources is discretized to create a grid of 29 × 12 points separated by 25 cm along x and y axes, leading to N = 348. Consequently, to compute the performance metric (25) the radius threshold deciding whether the localization is a true positive or not is xed to 25 cm.

To simulate a scenario indoor, the data is synthesized with the mirror source model, and a reection of the 3rd order. Then, to evaluate the eect of the model mismatch, A K is build by still considering a source-mirror model, with reections ranging from the 3-rd to the 0-th order.

The g. 4(a) plots results without model mismatch, i.e. with the 3-rd order of reection. All methods quickly reach 100% of true positives when the number of sensors increases. One can note that our approach is slightly better than the others, and SelRelax is the least good.

When the reection order decreases (g. 4(b), 4(c) and 4(d)), the model mismatch increases naturally, the performances degrade accordingly. One can note that our approach appears to be the less impacted, whereas SelGen and SelRelax seem to be the most sensitive to the mismatch. Also note that if the model error is too important (0-th reection order, i.e. free eld model), all methods performance are heavily impacted.

To complete the analysis via the True positive rate, the g. 5 plots the Jaccard index. It focuses on the proposed method and the random selection, as these are the two top performing methods according to g. 4. Both indicators have relatively similar trends, although the Jaccard index is slightly lower than the True positive with 0-th and 1-st reection orders in the model. In these two cases the dierence reveals the presence of False negatives and False positives.

From these studies, we may stipulate that our approach is more robust to model errors because it is the only data-informed method.

B. Results from experiments

For sake of comparison, we proceed with the same performance analysis as in previous simulations: methods are compared as a function of the number of sensors K, and the order of reections. The N r = 300 realizations for the Monte Carlo study are obtained based on the experimental setup described previously in sec. IV A.

The main dierence with sec.VI A is that we do not know how many reections are relevant to accurately model the received signal. One way to evaluate the order of reections is to compute the correlation coecient between the experimental and the synthesized data, based on the same conguration. The table I provides the average of theses coecients over the 25 available positions of the source, as a function of the reection order. As expected, increasing the order in the synthesized data leads to a higher correlation coecient, meaning that the mirror-source model is more accurate. This correlation seems to be maximized at the 2-nd order, and slightly decays back beyond.

Performance analysis

Based on the results in table I and to t with the graphs of sec. VI A, we display in the following the results with A K built for orders from 3 to 0. The plots are given in g. 6.

One can note that the results considering 3 or 2 reections are similar, which is in line with the correlation coecient table I. We suggest that this similarity could occur because the mirrored sources beyond the 2-nd order do not contribute to improve the source localization as we could expect. It can also be thought that other model mismatches (such as wall positions and orientations, exact source positions, etc.) aect the source-mirror model beyond a given order. Surprisingly, regarding of all methods, the performance trends reveal that the studied mismatch has a lower impact in the present real scenario.

Finally, apart from the model with no reection, our proposed method leads to better results than the ones from the state-of-the-art. Again, the Jaccard index is plotted in the g. First, the stopping criterion of the algorithm 1 is a chosen sensor number K. In practice, this remains a limitation because K may be too small to reach the asymptotic performance obtained when K = M . It also may be chosen too large, so that the addition of the latest sensors

1 :

 1 Output: K with |K| = K Initialize : take p k , p l ∈ M randomly, and initializeK = {p k , p l } ; while |K| < K doEstimate ΣK via eqs (11) and (12) ; Find p = arg max p∈M\K 1 + σ -2 a(p)ΣKa(p) H ; K ← K ∪ p; end Algorithm Greedy algorithm for data-informed and online sensor selection.

  The selection step requires an extensive search on M -|K| candidates: the computation cost slightly reduces as K is growing in size. An upper bound in complexity of this step is O(M N 2 ), and O(KM N 2 ) for the selection of K sensors accordingly.IV. CASE OF STUDY: INDOOR SOURCE LOCALIZATIONThis paper illustrates the method of sensor selection with an application to the near-eld acoustic source localization indoor. Experiments in situ are performed in a room equipped with a very large array. All the numerical and experimental results below follow the same common conguration, described in the next section. A. Experimental scenario: a 1020-microphone array mounted on walls The room and the antenna hardware are fully described in 10,11 . The picture in g. 1 shows the global setup during an acquisition. The array is made up of M = 1020 digital MEMS microphones, and is supported by the walls of the room of shape 8.01m × 3.75m × 2.94m. The sensor positions p m on these walls are random, as illustrated in g. 2(a). To ensure a tractable source localization, the source eld domain is restricted to the horizontal plane highlighted by the red laser in g. 1. The height of this plane is 1.35 m, and the acoustic sources (cf. the blue speaker, upper right) are positioned in the plane thanks to the laser level. An acquisition set is created at a sampling frequency of 25 kHz, for 25 source positions plotted with red points in g. 2(b). Each acquisition contains one single source in emission. In this way, it is easy to create exible source eld congurations by summing each of these acquisitions at will. The underlying motivation is to enable to do Monte Carlo simulations on real data. The source is assumed to be omni-directional and pointwise. At each acquisition, the speaker emits a sinusoidal component at the frequency f = 500 Hz. The xand y-axis true source positions are measured in reference to two vertical walls with as laser tape. To turn the source localization in the harmonic form, discrete Fourier transforms are made on frames of 6250 time samples (i.e. 0.25 s), to extract the chosen frequency from them. B. Propagation model indoor and model mismatch As a starting point, let us recall how to write the propagation model in a near and free elds:

3 .

 3 The number of mirrored points depend on the chosen order of reection. For instance, a model with reections: of the rst order needs the purple duplicates in g. 3(b); of the second order needs both purple and red duplicates.

FIG. 1 .

 1 FIG. 1. Experimental setup indoor: 1020 microphones mounted on walls. (color online)

FIG. 3 .

 3 FIG. 3. (a) Schema of the image-source model used to take into account the reections. The blue dot corresponds to the source position, the magenta dots are the rst-order positions, and the red dot is the second-order position. (b)The dierent discretized spaces in 2 dimension created up to the second order of reections. The blue space is the initial room, the magenta space is the rst order and red the second order. A 50cm gap was preserved between walls and the source grid. (color online)

  the article, Chepuri et al illustrate their studies with a source localization scenario based on time of arrival measurements. Since we use a dierent model, two key adjustments should be emphasized for a proper implementation. First, the problem must be overdetermined: the number of unknowns must be smaller than the number of sensors K. As such the linearized model (2) is not convenient when the source grid size N is large. As explained in 7 , it is tackled by staying with the continuous nonlinear form of the propagation model, and calculating the analytic FIM with 5 scalar parameters: the x-and ycoordinates of the source position, the source amplitude and phase, and σ 2 . It leads to a FIM of size 5 × 5. This choice of variables inherently describes the presence of only one source in the eld. Therefore, as the second key point, Chepuri et al develop an idea to discretize the domain of values potentially taken by the 5 parameters, while preserving the same overdetermined form. To do so, the relaxation problem replaces the use of a FIM by the use of a block diagonal matrix; each block is one FIM computed with a set of parameter values 7 (Sec. III.a). We use CVX toolbox 19 to solve the relaxed problem. To keep the computation tractable, the parameter domain discretization should be coarse because the size of the block diagonal matrix quickly grows. Among the 5 scalar parameters, we only discretize the x-and y-coordinate parameters at the 25 true source positions. It leads to a block diagonal matrix with 5 × 25 = 125 rows and columns.

  Correlation coecient between the experimental and the synthesized data, averaged on the 25 available source positions.

FIG. 4 .FIG. 6 .

 46 FIG. 4. Average percentage of true positives according to the number of sensors for our approach (blue line), SelRand (red line), SelGen (magenta plus) and SelRelax (black crosses). The data is synthetized with a 3-rd reection order and AK is dened for reection of the (a) 3-rd, (b) 2-nd, (c) 1-st and (d) 0-th order. (color online)

2 .

 2 7, with a focus on the iterative selections only. Similarly to the simulation results, the Jaccard index is similar to the True positive rate with the 2-nd and 3-rd reection order models. The dierence is marginal with 1-st order, but becomes signicant with the free eld model. Again, this dierence indicates a presence of False positives or False negatives that degrades the eective performance of source localization. In addition to the previous Monte Carlo analysis, the g. 8 provides an illustration of source localization results in one specic realization. It shows the absolute values of source vectors sK , estimated with the dier- ent methods, with S = 4 sources, K = 50 sensors, and a model with the 3-rd order of reection. The red dots correspond to the true source positions, the black `+' (resp. `x') signs provide the estimated source locations when they are true positives (resp. false positives). These gures reveal a sparse structure as expected, even if they do not contain 4 non-zero peaks only (then it can eventually leads to false positives). Our method (g. 8(a)) appears to be the most sparse. Although theses maps are specic to one realization, we veried that theses observations remain true for any other one. Statistical analysis of the selected sensors This last section investigates if there are recurring patterns of selection in the dierent Monte Carlo realizations. In other words, one may wonder whether some sensors of the array are more important than others, regarding the information they provide. Indeed, the proposed selection is data-informed and depends on each specic realization, i.e. the source positions. Nevertheless, does the selection prioritize some sensors more often in average? This question naturally leads to analyze the probability that the k-th sensor is selected among the r rst ones. As such, it relates with the following empirical cumulative distribution (cdf ): P(p k ∈ K; |K| ≤ r).

( 27 )

 27 Obviously, this probability monotonically increases since it is a cdf, and is dened for 2 ≤ r ≤ M . The g. 9 draws this empirical cdf for each of the M = 1020 sensors. It evidences that cdfs are heterogeneous, and reach 1 at a dierent rate according to the sensor index k. Consequently, some sensors are more likely to be selected than others to create an array of r (or less) elements. To illustrate this, we display in g. 10, the 60 sensors having the highest values P(p k ∈ K; |K| ≤ r) for r = 50, and the color indicates the corresponding value.Interestingly, most of these sensors lie in the plane in which the sources belong to.RemarkAt this point, the statistical study indeed reveals recurring selection patterns in the dierent realizations. One intuition would be to extract the most appearing sensors to design one global array. In this way, array design could be done once and oine for a realistic scenario, thanks to the proposed online sensor selection performed on a set of simulated data. Although beyond the scope of this paper, it seems to be an interesting future investigation. VII. DISCUSSION AND OPEN PROBLEMS With this paper, the authors wish to initiate the investigation of source localization helped by iterative sensor selection. The presented work experimentally validates the interest of this research direction, however some problems remain open.

  brings a marginal improvement on the source location estimation. Thus, an alternative stopping criterion based on an automatic choice of K is a future direction to investigate in order to enhance the use of the algorithm.The second open problem deals with the computation time. Indeed, SBL needs to be repeated K times, and has a complexity of order O(|K| • N 2 ) 45 . If K ≈ M it is trivial that the iterative approach requires more computation than doing SBL once with the full array. Nevertheless, the computation time of the algorithm 1 is generally more time-consuming, since the dimension of the source vector N mainly drives the complexity of SBL. Replacing the conventional SBL with an accelerated technique of lower complexity is consequently one interesting direction to explore see for example the state of the art introduction in29 .

FIG. 8 .FIG. 9 .

 89 FIG. 8. Source localization estimation done by SBL sK for 1 realization, with S = 4 sources and K = 50 sensors. (a) Proposed approach, (b) SelRand, (c) SelGen, (d) SelRelax. (color online)

J. Acoust. Soc. Am. / 19 September 2022 Sequential sensor selection for sparse source localization

J. Acoust. Soc. Am. / 19 September 2022 Sequential sensor selection for sparse source localization


Christof Weber
email: christof.weber@phlu.ch

Janka Medova
email: jmedova@ukf.sk

Maryna Rafalska
email: maryna.rafalska@univ-cotedazur.fr

Ulrich Kortenkamp
email: ulrich.kortenkamp@uni-potsdam.de

Simon Modeste
email: simon.modeste@umontpellier.fr

Introduction to the papers and posters of TWG11: Algorithmics

Keywords: Algorithms, algorithmic thinking, teaching and learning of algorithms

In CERME12, our working group "Algorithmics" started its work as a newly established TWG. Since algorithms have always been at the heart of mathematics and their importance has been steadily increasing since the beginnings of theoretical computer science, the design and analysis of algorithms -called algorithmics [START_REF] Traub | Iterative methods for the solution of equations[END_REF][START_REF] Knuth | Algorithmic thinking and mathematical thinking[END_REF]) -lies at the intersection of mathematics and computer science. For this reason, on the one hand, various algorithms and algorithmic activities have their traditional place in mathematics curricula at all levels. At the school level, mathematics and computer science have interacted since the 1980s, when many schools set up labs with computers equipped with programming software. On the other hand, many questions arise in the context of teaching and learning algorithms: a first, more applied group of questions aims at algorithms in mathematics education and curricula, a second, more theoretical group of questions seeks to clarify the concepts of algorithm and algorithmic thinking.

Conference presentations

Due to the Corona pandemic, the conference was held as a virtual event. Nevertheless, a total of 11 papers and 7 posters were presented remotely by their authors at the conference, with a total of 24 group participants from 11 countries. The contributions were considered in four themes, as follows.

Theme 1: Beliefs and domains in which algorithmic thinking occurs

A first group of papers focuses on the place and importance of algorithms in mathematics in general and arithmetic in particular. They assess the beliefs of experts about the role of algorithms in mathematics and mathematics education or their role in mathematics courses.

• Lockwood, DeJarnette, Thomas and Mørken offer three perspectives on algorithms, particularly in computational settings: an algorithmic approach in a mathematical example, the view of a mathematician, and the view of an undergraduate students taking a course in mathematics. • Geraniou and Hodgen interviewed two mathematics educators who had experience using technology to solve mathematical problems, and they, too, shared very different views on algorithms in mathematics education, one even not seeing the use of algorithms as a mathematical activity. • Kortenkamp analyzes an arithmetic course for pre-service primary teachers. He identifies several algorithmic activities in the topics covered in the course, such as designing algorithms, specifying algorithms, performing algorithms, proving their correctness, and comparing algorithms.

• Leifeld and Rezat's poster provides a thorough analysis of the possibilities of certain arithmetic algorithms for addition and subtraction to deepen students' understanding of inverse operations.

Theme 2: Teaching and learning of algorithmic thinking at primary level Another group of papers focuses on teaching and learning algorithmic thinking in primary school.

They use different tasks with different goals: Some use algorithmic thinking as a means to an end (in the sense of learning a new mathematical concept), some use algorithmic thinking as a goal (in the sense of understanding a given algorithm, or developing an algorithm to solve a problem):

• Crisci, dello Iacono and Ferrara Dentice explore how primary school children can be stimulated to learn new mathematical concepts by Scratch. In their report, they present a specially designed task that required visual programming to complete a given figure so that it becomes axially symmetric. The children developed two different strategies to solve this task. For example, they found out that points that are axially symmetrical to each other must be equidistant from the axis. • Funghi and Ramploud are interested in how to teach the standard long-division algorithm so that children understand why it "works." To this end, they had fourth graders compare the optimized, digit-by-digit long-division procedure with the procedure in which the divisor is repeatedly subtracted from the dividend. Their analysis of class discussions suggests that this approach could actually result in less rote learning, but in more conceptual learning. • Zindel's study wants children to acquire algorithmic thinking without using computers. In her papers, children are instructed to decrypt and encrypt certain words. They had to articulate the necessary steps themselves and record them in writing. Although these texts show great differences, the author succeeds in reconstructing some constituents of algorithmic thinking. • In Gaio's study, too, children are asked to develop algorithms in the sense of systematic procedures, without the help of computers. Here, the tasks given to children of different school levels (3rd to 8th grade) are in the context of sorting problems. As the author reports, he can see traces of classical algorithms in the procedures that the children have worked out cooperatively.

Theme 3: Teaching and learning of algorithmic thinking at university level

Concerning the development of algorithmic thinking at university, the discussions showed two big issues: the development of algorithmic thinking and algorithmics in mathematics, for students, independently of their projects, and more specifically, the development of algorithmic thinking in mathematics for future teachers, and in particular future primary teachers.

Four papers dealt with algorithmic thinking in advanced mathematics, three at university, and one concerning an education program for gifted students. Above them, three were interested in links with discrete mathematics, combinatorics, graph theory, which illustrates the specificity of those mathematical fields, at the interface with computer science:

• De Chenne and Lockwood explore the use of programming and computer science in solving basic counting/combinatorics tasks in college, and how the knowledge of student in computer science can influence their solving strategies and support their learning. • Medová, Milicic and Ludwig study the competencies involved in algorithmic thinking for university, and in particular abstraction, modelling, and visualization skills which are difficult to master for students, and questions the development of computational thinking in mathematics.

• Bóra and Gosztonyi analyze the place given to algorithms and what could be seen as algorithmic problem solving, in Hungary's advanced mathematics programs, questioning what can be considered traditionally as algorithmic in mathematics and its place according to mathematical culture of the country. • The paper of Calor, Palha and Kubbe concerns at the same time advanced mathematics and preservice secondary teachers' education. It deals with analysis, and in particular developing instructional material concerning differential equations for algorithmic thinking and programming. First results show that students did indeed develop algorithmic thinking in their work.

The two other contributions dealt with algorithmic thinking in and for primary teacher training: • Weber's paper examines primary teachers' use of loops to solve a geometrical problem and their conceptions of the loop construct. It elaborates some challenges in their conceptions and some misconceptions that require deepening their understanding from a teacher training perspective. • Dobgenski and da Fontoura's poster presents and reflects on an experience of making pre-service primary teachers deal with computational thinking using Scratch.

Theme 4: Concepts related to algorithmic thinking: computational thinking, algebraic thinking, problem solving, and mathematical literacy

The last group of contributions deals with different, no less relevant aspects of algorithmic thinking: • Rafalska's paper illustrates how tasks could be constructed in order to lead children in mathematics lessons (without the use of computers) to algorithmic thinking in the sense of developing a solution strategy and which individual learning processes can be triggered by these tasks. • Pohlkamp and Lengnink's paper takes a different look at algorithms: It discusses algorithms that make decisions and are thus socially relevant. Addressing and studying them in the classroom would mean taking more seriously the educational mandate to teach social skills as well.

Finally, two poster proposals deal with two concepts related to algorithmic thinking:

• Rekstad and Rasmussen investigate the question to what extent aspects of computational thinking mentioned in the literature are also reflected in teachers' beliefs when asked about the role of computational thinking in mathematics education. • The relationship between algorithmic and algebraic thinking is the subject of Müller-Späth, who plans to investigate how algorithmic thinking (realized by an app) affects the development of the ability to generalize and thus of algebraic thinking.

Conference discussions

As mentioned earlier, our working group has just begun its work, and a common understanding of the concepts has yet to be developed: What does algorithmics mean in the context of teaching and learning mathematics? What is algorithmic thinking? To this end, after the presentations in which quite different views were expressed, we worked on the following three questions:

Question 1: Which mathematical algorithms could stimulate algorithmic thinking?

The discussion of this question revealed relatively unanimously five mathematical types of algorithms: i. Algorithms based on the place value system (standard algorithms for addition etc., algorithm for calculating logarithms), ii. graph-theoretic algorithms (shortest path problem, Königsberg problem), iii. approximation algorithms (Heron's algorithm, Newton's method, etc.), iv. sorting algorithms (heap sort, bubble sort, etc.), and v. miscellaneous (Gauss's Easter algorithm, etc.). We were not in agreement of whether each procedure is also an algorithm. For example, everyday procedures (tying shoes, making jam sandwiches, etc.) were not considered by all participants to be suitable for addressing and promote algorithmic thinking in its "proper sense" because they show only one characteristic feature of algorithms: the order of steps.

Question 2: Which mathematical topics could promote algorithmic thinking?

The discussion of this somewhat broader question also yielded five topics from which tasks could come to stimulate algorithmic thinking: i. number theory (arithmetic, prime number tests, factorization), ii. discrete mathematics (graph theory, combinatorics, counting problems, etc.), iii. geometry (transformations, algebraic geometry), iv. computer science (cryptography, etc.), and v. games and puzzles (Rubik's cube, tower of Hanoi, etc.). One participant's question about what properties these fields would have in common was discussed intensively and controversially.

Question 3: What (human) activities with algorithms can we think of?

The activities discussed suggest a wide range of possible activities to deal with algorithms: i. creating (developing algorithms, improving algorithms, debugging algorithms, optimizing algorithms, transferring algorithms to an analogous situation, etc.), ii. analyzing (effectiveness and proof, efficiency, complexity, stability, similarity etc.), and iii. comparing (comparing different algorithms for the same problem, comparing analogous algorithms for different problems, classifying algorithms, etc.). Although executing an algorithm without any reflection would be a possible activity with algorithms, most participants do not want this to be understood as algorithmic thinking.

Surely the reader can think of further examples or answers to these questions. In other words, the three questions need to be discussed further and their answers are still quite open.

Outlook

As the overview of the contributions as well as the first answers to central questions show, there is a great variety of approaches (theories, methods) and views (topics, perspectives) in our working group. Given that we are entering a young (or at least long-neglected) area of research in mathematics education and that we have just begun work in our TWG, it was to be expected that the results would be disparate and sometimes controversial. However, in terms of a first step towards a robust and sustainable understanding of concepts, this diversity makes us confident that there are many more questions around the challenging topic of algorithms and algorithmic thinking that are worth working on. With this in mind, we look forward to CERME13 and hope for a fruitful continuation of the work we have begun -and that it can then be carried out again as a physical conference.


Renaud Chorlay
email: renaud.chorlay@inspe-paris.fr

Tanja Hamann
email: hamann@imai.uni-hildesheim.de

Jenneke Krüger

Antonio Oller-Marcén

CERME 12: Thematic Working Group 12 History in Mathematics Education

Introduction

History of mathematics in mathematics education, and history of mathematics education continue to receive much attention. Although empirical research and coherent theoretical/conceptual frameworks within this area have emerged relatively recently, there exists an increasing interest in these lines of work. The purpose of this CERME TWG is to provide a forum to approach mathematics education in connection with history and epistemology.

TWG12 welcomes both empirical and theoretical research papers, and posters proposals related to one or more of the following issues -although any paper or poster of relevance to the overall focus of the group is taken into consideration. All educational levels can be considered, from early-age mathematics to tertiary education and teacher training:

1.
Design of teaching/learning materials using the history of mathematics, preferably with conclusions based on empirical data.

2.

Research on the existing uses of history or epistemology in curricula, textbooks, or classroom practice.

3.

Research on the history of mathematics education, on local, national or international level.

4.

Connections between mathematics education, history of mathematics, and history of mathematics education: Theoretical and methodological issues.

Even though the inception of this TWG is fairly recent -it started in CERME6 ( 2009) -it has deeper institutional roots within the mathematics education research community. The History and Pedagogy of Mathematics (HPM) study-group (http://www.clab.edc.uoc.gr/hpm/about%20HPM.htm) was created at the 1972 ICME conference; it has been organizing satellite conferences to the ICME meetings since 1984, and has several active regional branches (HPM-Americas, European Summer Universities).

In CERME12, TWG12 accepted 13 papers and 1 poster. One paper was withdrawn due to financial problems, one participant could not present because of illness. The Covid pandemic was the reason that CERME12 was delayed until 2022, however, even in 2022 all sessions were on-line. Thanks to the excellent virtual platform created by the Free University of Bolzano and the effort and enthusiasm of the organization and participants, TWG12 was again a success, though with participants from a more restricted range of countries than in previous years Before going into any details, it should be stressed that TWG12 has four general but distinctive features which give these meetings their specific flavour. Firstly, its topic lies at the intersection of different fields of research -mathematics education research and history of mathematics -which requires versatility and methodological vigilance [START_REF] Fried | Can Mathematics Education and History of Mathematics Coexist[END_REF][START_REF] Chorlay | History of Science, Epistemology and Mathematics Education Research[END_REF]. Secondly, the strength of the historical and the HPM community varies greatly among countries, and these meetings play a crucial role for researchers working in relative isolation, and with difficult access to resources in the field. Thirdly, the scope of TWG covers both history in mathematics education and history of mathematics education, which are two significantly different research topics (TSG 27 and 55 in ICME14); connecting the two lines of investigations is a constant challenge. Fourthly, since the topic of TWG12 is neither specific to one level of the educational system (from primary education to teacher-training) nor to any single mathematical topic (be it fraction concepts, algebra, proof, etc.), the work in TWG12 intersects that of most other TWGs. It should be noted that, as in former editions of CERME, there was little intersection with what was covered in TWG8 (Affects and the teaching and learning of mathematics), TWG10 (Diversity and mathematics education), in spite of the fact that it is not uncommon for outsiders of the HPM research community -among which most policy-makers and curriculum-designers -to ascribe such goals to the historical perspective in teaching.

These four features make these meetings not only useful but also challenging and exciting. As the final discussion made clear, the general feeling among the participants was that one of the main outcomes of this meeting is that we actually learned a lot from the one another, both from their papers and from the lively discussions.

Some significant features of the 2022 conference

From an organizational viewpoint, due to the relatively small number of contributions, the TWG12 team decided to offer three work formats for the sessions. In addition to the two traditional formatsoral presentations of a paper (talk) and discussions of topics of general interest -we also organized two 1h30-hour workshop sessions (with two parallel presentations each). A distinctive feature of TWG12 is that it potentially covers all mathematical topics, all historical periods, and many different research perspectives (history of education, historical roots of a given concept or practice, analysis of teaching resources and curricula, task design, analysis of teaching practices, training design etc.), usually making it difficult for the participants to make the most of short oral presentations and 8-page papers. The workshop format enabled some of the participants to actually share some material, and discuss both research questions and methodological issues. Out of the 13 contributions four were found particularly suitable to be presented in a workshop: two contributions on history of mathematics education (Hamann & Schmidt-Thieme, Goemans & De Bock) and two on history in mathematics education, either at the secondary level (Spies & Junker) or in teacher education (Arnal-Bailera, Beltrán-Pellicer & Oller-Marcén). The participants of TWG12 agreed that the new format was enjoyable and should be kept for upcoming CERME conferences, whenever possible.

As to content, there was a general balance between the two main lines of investigation: 6 papers on the history of mathematics education, and 7 papers on history in mathematics education. Let us highlight a few features of the Bolzano edition of CERME, without aiming for a thorough survey.

As far as history of mathematics education is concerned, a significant share of the papers bore on the New Maths movement in Europe, thereby shedding light on this complex historical phenomenon from a variety of perspectives: connections with advanced mathematics (in terms of content and of interactions among different social groups and stakeholders); controversies at the time of the reform and long-term effects (if any); impact on curricula, teaching resources, and actual teaching practices; specifics of different countries (Belgium, Nordic countries) and of different segments of the education system (e.g., vocational training).

As far as history in mathematics education is concerned, the "traditional" line of work on task design was represented by several papers (Surroca, Spies, Oswald). In spite of differences in terms of mathematical contents and segments of the educational system at stake (from primary to tertiary education), these papers illustrate how the HPM perspectives requires that researchers combine inputs from two different fields of research, mathematics education and history of mathematics. They also showcase a practice which lies at the heart of the HPM perspective: the use of original historical sources in the classroom [START_REF] Chorlay | Historical sources in the classroom and their educational effects[END_REF].

Over the last decade, teachers have become a central topic of interest for researchers working in the HPM community, with inputs from many theoretical frameworks. Teacher education was the main focus of the papers presented in CERME 12, with analyses making use of a variety of theoretical constructs and frameworks from mathematics education: the "four S"-model (van den Bogaart), Zhang's the ADTRE model, and Godino's "epistemic configurations" (in the paper of Oller-Marcén et al.). By contrast, the analysis of "ordinary" teaching practices in primary or secondary education is still an emerging topic, with one paper in CERME 12 (Chorlay), with theoretical inputs from the Documentational Approach to Didactics (Gueudet, Pepin, Trouche) and the Didactic and Ergonomic Double Approach (Robert, Rogalski, Vandebrouck).

Challenges for the future

The number of papers and posters presented in TWG12 was modest, and half of the participants of this edition of TWG12 were first-timers. This reflects the fact that, in spite of the fact that the HPM community has a stable international basis and ramifications in many countries, it remains somewhat marginal within the mathematics education community. Some of the reasons accounting for this situation are structural and beyond the scope of our actions, such as the modest role of history of mathematics in curricula, and the necessity for researchers in this field to be knowledgeable in two different areas (history and mathematics education). However, the importance of the inclusion of history of mathematics in teacher training seems to be gaining more recognition in many countries, thus providing opportunities to further and disseminate research. Also, the fact that research carried out from an HPM perspective touches on virtually all mathematics education topics should pave the way for exchanges and collaboration. While this is true at a theoretical level, it remains a challenge to strengthen and stabilize these connections. In the spirit of CERME's three Cs (Communication -Cooperation -Collaboration), it is one of the aims of TWG12 to reach out and generate fruitful interactions.