F Plourabou
email: fplourab@imft.fr

P Uszes

R Guibert

F Plouraboué

Source identification of propagating waves inside a network

 to an arbitrary number of sensors, and estimate its complexity which depends on sensors number Ns and time dicretisation Nt. Increasing detectors number increases source identification robustness to noise up to a cut-off number being a small fraction of N . This cutt-off detector density for efficiency in noise reduction is of practical significance. We also analyze and discuss the method sensitivity to total recording time Te, sampling frequency, and signal to noise ratio. Finally we propose a pre-sectorisation to improve the systematic exploration algorithm proposed in [1] and we show a O(N ) drop in complexity leading to a O(N 3 ) cost for the source identification. Several tests of the method on model and real graphs confirm and support the presented results.

INTRODUCTION

General context and motivations

In various contexts like water-hammer waves within pipes networks, river tidal waves in river basin, explosives acoustic waves in gas pipes tubing, electrical waves in powerdistribution networks, it is interesting to locate the origin of waves propagating inside networks [START_REF] Caputo | Inverse source problem in a forced network[END_REF], [START_REF] Caputo | Inverse source problem in a forced network[END_REF], [START_REF]Traffic Flow on Networks[END_REF], [START_REF] Bressan | Flows on networks: Recent results and perspectives[END_REF], [START_REF] Giudicianni | Topological Taxonomy of Water Distribution Networks[END_REF], [START_REF] Prakash | Review of power system distribution network architecture[END_REF], [START_REF] Wishart | Advances in metabolite identification[END_REF], [START_REF] Caputo | Oscillations of networks: The role of soft nodes[END_REF]. The simplest exemple of spatio-temporal source identification is the localization of a sudden breakage or leakage event within a single pipe [START_REF] Gao | A model of the correlation function of leak noise in buried plastic pipes[END_REF]. In this case, the origin of wave event propagating at velocity c, in both ±x directions, is searched from two distinct sensors recording the signal at some location along the pipe axis x. Standard signal processing methods can then be applied [START_REF] Gao | A model of the correlation function of leak noise in buried plastic pipes[END_REF], looking for the highest correlation between temporally translated signals at each sensor location within the two-dimensional half-space x -t. Nevertheless, even in the most simple networks, wave rebounds, or multiple path-waves produce as many temporal patterns for the signal recorded at the sensor location as possible wave-paths from event to sensor preventing simple signal processing methods to provide an answer to locate sources. Networks are hereby defined by their non-oriented graph G(V, E), vertex (node) set V , and edge set E, having a symmetric adjacency matrix. The wave propagation along each edge is supposed to arise at constant (possibly non-homogeneous) speed. The wave spreading into the network, i.e its bouncing and propagation at each node are determined by the graph Laplacian matrix ∆ G which sets the network's harmonic coupling between vibrating modes [START_REF]Traffic Flow on Networks[END_REF], [START_REF] Bressan | Flows on networks: Recent results and perspectives[END_REF]. For general networks four important questions are interesting to address :

Manuscript received 04-Feb-2021 (i) For a (non dissipative) wave equation is there a way to find the -possibly unambiguous-source location in time and space from recording the induced signal at some sensor location ? (ii) How is the identification method threatened by noise and possibly improved by increasing the number of sensors ? (iii) Is there some possibly optimal way of distributing the sensors inside the networks so as to improve both identification and its reliability ? (iv) Finally, is it possible to improve the computational efficiency of source identification ?

Here we consider 'single event' localization only as in [START_REF] Caputo | Inverse source problem in a forced network[END_REF], [START_REF] Shen | Locating the source of diffusion in complex networks by time-reversal backward spreading[END_REF], [START_REF] Wang | An algorithm for locating propagation source in complex networks[END_REF], [11]. Question (i) has been addressed in [START_REF] Caputo | Inverse source problem in a forced network[END_REF] for wave-source identification in the case of two sensors. We provide an extended description of their approach in §2.1. Furthermore, in the context of diffusive processes [START_REF] Shen | Locating the source of diffusion in complex networks by time-reversal backward spreading[END_REF] found an identifiability condition associated with the need for signal not having the same time to reach any detector from two distinct sources. More general source location algorithms for heterogeneous propagation models have been proposed [START_REF] Wang | An algorithm for locating propagation source in complex networks[END_REF], [11].

Issue (ii) has been partly addressed in [START_REF] Caputo | Inverse source problem in a forced network[END_REF] in the case of two sensors for intern-type event identifications. Also in [START_REF] Shen | Locating the source of diffusion in complex networks by time-reversal backward spreading[END_REF] the influence of random time-delay has been investigated for some discrete value of noise. Question (iii) has been qualitatively addressed in some contributions, using networks concepts such as node's centrality [START_REF] Scardoni | Centralities Based Analysis of Complex Networks[END_REF], [START_REF] Yang | Locating the propagation source on complex networks with Propagation Centrality algorithm[END_REF], [START_REF] Hongjue | An universal algorithm for source location in complex networks[END_REF]. Concerning (iv), in the context of diffusive processes over graphs source identification complexity has received dedicated attention in the literature [START_REF] Shen | Locating the source of diffusion in complex networks by time-reversal backward spreading[END_REF], [START_REF] Brockmann | The hidden geometry of complex, network-driven contagion phenomena[END_REF], [START_REF] Peng | Smartphone malware and its propagation modeling: A survey[END_REF], [START_REF] Pinto | Locating the source of diffusion in large-scale networks[END_REF]. In the case of breadth-first tree structure of first arrival, [START_REF] Pinto | Locating the source of diffusion in large-scale networks[END_REF] have provided a O(N 3 ) cost method for source identification. Furthermore, in the case of diffusion propagation within a network [START_REF] Shen | Locating the source of diffusion in complex networks by time-reversal backward spreading[END_REF] have proposed a O(N 2 ln(N )) method using a time-reversal backward spreading algorithm. In some applications, the 'event' can be described by the presence of an external source term appearing at given nodes of the network (subsequently called an intern-type event), as in the case examined in [START_REF] Caputo | Inverse source problem in a forced network[END_REF]. Nevertheless, in other applications, the event results in a structural change of the network itself for which a given node becomes exclusively triggered by an external source impinging its interactions with neighbor nodes of the network. This is found in network pipe disruption, for which the external pressure is suddenly applied at one node, whatever the pressure nearby. Those boundary-condition driven events (subsequently called BC-events) produce a subtle change in the graph topology illustrated in Figure 1a-b. They suddenly locally shift the graph topology from locally un-directed (Cf Figure 1a) to locally directed (Cf Figure 1b). This subtle change has a mostly unfortunate consequence : the graph mathematical structure is changed by BC-events, and the results provided in [START_REF] Caputo | Inverse source problem in a forced network[END_REF] no more apply in this case. Hence, the aim of this contribution is to -Extend the identification algorithm proposed in [START_REF] Caputo | Inverse source problem in a forced network[END_REF] to BC-events identification with an arbitrary number of sensors. -Analyze and discuss the method sensitivity to sensor's recording frequency, total recording time T e and signal to noise ratio. -Propose a pre-sectorisation method reducing the numerical cost of source identification.

Section 1.2 provides some useful definitions. Section 2 provides the mathematical analysis of two distinct extensions of [START_REF] Caputo | Inverse source problem in a forced network[END_REF], the details of which is first presented in §2.1. Then, the source identification proof for an arbitrary number of sensors is detailed in §2.2. Lately, §2.3.1 provides the theory for BC-events identification. Section §3 describes the algorithms for both intern type or BC-type source identification. In §4 practical improvements are discussed, among which a pre-sectorisation method to improve the algorithmic cost is presented in 4.3. §5 provides numerical tests and illustrations of source identification in various real-world graphs.

A nomenclature is provided in Appendix A.

Some concepts and definitions

Definition 1 (Source).

A source is the location s ∈ V , within the vertex set of the graph, of a time-varying perturbation, an event, arising at a single node. The event can either be of intern type or of Boundary-Condition type (i.e BC-type). The type of the event also define the type of the source, i.e intern source or BC-source.

Definition 2 (Event).

An event is defined by the couple (λ(t), S), where λ(t) is an applied, compact support, temporal perturbation λ(t) ∈ L 2 and S a source vector, with S ∈ {0, 1} N , S = ŝ, i.e, having components S k = δ sk . The event can either be of intern type or of Boundary-Condition type (i.e BC-type).

Definition 3 (Intern event

). An event is said of intern type if its contribution is only added from the contribution of neighbors of this node. internal symmetries. From figure 2, it is then easy to find that some sources are not possible to distinguish with limited detector set. The fact that node #3 of figure 2 does not permit to built a strategic set can also be found from computing the eigenvector set of the graph Laplacian given in [START_REF] Caputo | Inverse source problem in a forced network[END_REF], where we showcase the third line (associated with node #3 highlighted in red) which violates Definition 6 for strategic set.

V =       0.4472 -0.6533 0.
      (1)

THEORETICAL ANALYSIS

This section details Caputo et al. [START_REF] Caputo | Inverse source problem in a forced network[END_REF] results in 2.1 the extension of which for an arbitrary number of sensors is detailed in 2.2. Furthermore, the extension of Theorem 1 for uniqueness of BC-events identification is given in 2.3.

Detailed presentation of [1]

2.1.1 Evaluation of X.

The wave amplitude at time t ∈ [0, T e ] at each vertex of graph G(V, E), with card(V ) = N , is represented by vector

X X(t) = (X 1 (t), X 2 (t), . . . , X N (t)) (2)
With T e the final recording time. ∆ G is the Laplacian matrix of non-oriented graph G. The governing equation of X(t) is given by the d'Alembert operator

Ẍ -∆ G X = λ(t)S for t ∈ [0, T e ] X(0) = a ∈ R N and Ẋ(0) = b ∈ R N (3)
where a et b are the initial state and velocity of the system. λ is the time-variation of a single source localized in single node (vertex) s ∈ V , with S i = δ is . Model (3) represents forced harmonic vibrations into a network resembling coupled vibrating strings, with a time-varying force imposed at node s.

It is assumed that (λ(t), S) are chosen in admissible event sets defined by

A s :=    ∃T 0 ∈ [0, T ]/λ(t) = 0 ∀t ∈ [T 0 , T ] S = (s 1 , s 2 , . . . , s N ) with s k ∈ {0, 1} k ∈ V, & N k=1 s k = 1 (4)
This admissible set provides the condition for vector source S to represent only one node among all possible ones, and the need for the forcing term to have compact support so that it disappears after finite time T 0 . The observation operator is defined in [START_REF] Caputo | Inverse source problem in a forced network[END_REF] as an application from the event space to observation space :

M (λ(t), S) = {D(t) with t ∈ [0, T]} (5)
where the detection vector D, is chosen at two sensors locations (i, j) ∈ V [START_REF] Caputo | Inverse source problem in a forced network[END_REF]. This abstract observation operator is useful to consider since it provides the formal link between the event, i.e the source node where some transient forcing occurs, and the later observation of the resulting disturbances within the network at some pre-defined locations provided by D. It is assumed that the eigenvalue of ∆ G are distinct and ordered

0 = -ω 1 2 > -ω 2 2 > • • • > -ω N 2 (6)
Since ∆ G is a symmetric matrix of a non-oriented graph, its eigenvectors {v 1 , v 2 , .., v N } produce an orthogonal base of R N , so that

X(t) = N n=1 y n (t)v n with y n (t) = X(t), v n ∀n ∈ [1, N ] (7
) Lemma 1. For T ∈ (T 0 , T ), and X solution of

Ẍ(t) -∆ G X(t) = 0 ∀t ∈ (T , T ) X(T ) ∈ R N et Ẋ(T ) ∈ R N (8)
ξ = (k 1 , k 2 , . . . , k Ns ) a strategic set. Then,

x k1 (t) = • • • = x k Ns (t) = 0, ∀t ∈ (T , T ) ⇒ X(T ) = Ẋ(T ) = 0 (9)
This lemma states that, if the restriction of the wave amplitude at strategic set (Cf definition 6 and figure 2 for illustration) location X|ξ(t) is zero at any time after the event to occur, then, the wave amplitude and velocity are zero everywhere at final time. Using lemma 1 [START_REF] Caputo | Inverse source problem in a forced network[END_REF] have shown that Theorem 1. If (λ, S) is part of admissible events set, i.e (λ, S) ∈ A s , X(t) solution of (3) given at nodes of strategic set ξ, i.e X| ξ (t), permit to find a unique X(T ).

The symmetric Laplacian matrix ∆ G admit an orthogonal base of eigenvectors. Considering the linear homogeneous differential problem Ẍ -∆ G X = 0, it is possible to find the general solution projected in the ∆ G eigenvector base V

X(t) = (y 1 (T ) + (t -T ) ẏ1 (T )) v 1 + N n=2 y n (T ) cos(ω n (t -T )) + ẏn (T ) ω n sin(ω n (t -T )) v n (10)
This projection permits to infer the signal everywhere in the network, from recording it in a limited number of sensors located in ξ ⊂ V . Components y n (T ) for n ∈ [1, N ] permits to built the mode amplitude vector y(T ) ∈ R N which gives, for each mode, the oscillation amplitude. Sensors at positions ξ indeed permit to evaluate the mode amplitude vector y(T ) and its time derivative ẏ(T ), encapsulated into 2N-vector Y(T ) ≡ [y(T ), ẏ(T )]as in [START_REF] Plouraboué | Identification of internal properties of fibres and micro-swimmers[END_REF]-from linear regression with observations. This linear regression uses the detection vector

D ≡ X k (t), k ∈ ξ with card(ξ) = N s , sampled at discrete value t m = m∆t, with m = 0, . . . , M D = X k1 (T 0 ), . . . X k1 (t m ), . . . , X k1 (t M ), X k2 (T 0 ), . . . X k2 (t m ), . . . , X k2 (t M ), X k Ns (T 0 ), . . . X k Ns (t m ), . . . , X k Ns (t M ) (11)
The regression aims minimizing X|ξ -D 2 , the difference between detection and the wave amplitude prediction at sensor location ξ, which is a function of unknown wave amplitudes Y(t). Since X fulfills a linear dependence with Y(t), a matrix R can be built such that, X|ξ = RY with regression matrix R equals to

              v 1 1 (T 0 -T )v 1 1 cos (ω 2 (T 0 -T ))v 2 1 sin (ω 2 (T 0 -T )) ω 2 v 2 1 . . . cos ω N (T 0 -T ) v N 1 sin ω N (T 0 -T ) ω N v N 1 . . . . . . . . . . . . . . . . . . . . . v 1 1 (t M -T )v 1 1 cos ω 2 (t M -T ) v 2 1 sin ω 2 (t M -T ) ω N v 2 1 . . . cos ω N (t M -T ) v N 1 sin ω N (t M -T ) ω N v N 1 v 1 2 (T 0 -T )v 1 2 cos (ω 2 (T 0 -T ))v 2 2 sin (ω 2 (T 0 -T )) ω 2 v 2 2 . . . cos ω N (T 0 -T ) v N 2 sin ω N (T 0 -T ) ω N v N 2 . . . . . . . . . . . . . . . . . . . . . v 1 2 (t M -T )v 1 2 cos ω 2 (t M -T ) v 2 2 sin ω 2 (t M -T ) ω N v 2 2 . . . cos ω N (t M -T ) v N 2 sin ω N (t M -T ) ω N v N 2 . . . . . . . . . . . . . . . . . . . . . v 1 Ns (T 0 -T )v 1 k Ns cos (ω 2 (T 0 -T ))v 2 Ns sin (ω 2 (T 0 -T )) ω 2 v 2 Ns . . . cos ω N (t 0 -T ) v N Ns sin ω N (t 0 -T ) ω N v N Ns . . . . . . . . . . . . . . . . . . . . . v 1 Ns (t M -T )v 1 Ns cos ω 2 (t M -T ) v 2 Ns sin ω 2 (t M -T ) ω N v 2 Ns . . . cos ω N (t M -T ) v N Ns sin ω N (t M -T ) ω N v N Ns               (12
) Then, the linear regression permits to find unknown wave amplitudes Y from

Y = [R • R] -1 R D (13)
and then X(t) is found from using ( 13) in [START_REF] Wang | An algorithm for locating propagation source in complex networks[END_REF].

Source location detection:

in order to find the source location, the first step consists in simplifying the temporal ordinary differential equation from projecting it into a suitable base. Since this temporal problem fulfills a Sturm-Liouville problem

φm (t) = µ m ϕ m (t) ∀t ∈ [0, T ] ϕ m (0) = ϕ(T ) = 0 (14)
with eigenfunctions ϕ m associated with eigenmodes µ m

ϕ m (t) = 2 T sin mπ T t et µ m = mπ T 2 (15)
Functions {ϕ m } form a complete orthogonal set of L 2 (0, T ) equiped with scalar product f, g = T 0 f (t)g(t)dt. This base is then used to simplify the temporal o.d.e problems of the wave amplitude (3). Multipling (3) by ϕ m and integrating twice by part leads to

-(∆ G + µ m I)X m = λ m S + P m , (16)
Where X m ≡ X, ϕ m , resp. λ m ≡ λ, ϕ m , are the projection of X, resp. λ, over ϕ m , and P m = φm (T )X(T )φm (0)X(0). Defining matrix

A m ≡ ∆ G -µ m I, (17)
(16) reads -A m X m = λ m S + P m . (18)
(18) thus represents the second step of the source location strategy for providing a necessary linear algebraic problem for projected mode m. In this formulation of the amplitude projection dynamic, the right-hand-side term P m takes care of the influence of initial and final conditions of both basefunction velocity φm and amplitude vector X. Then [1] search for [START_REF] Plouraboué | Identification of internal properties of fibres and micro-swimmers[END_REF] at every nodes expect those where sensors are located. Denoting (i, j) the set of the two sensors ( [START_REF] Caputo | Inverse source problem in a forced network[END_REF] consider the case of two sensors only), and removing the corresponding components on the right-hand-side leads to

-A (i,j) m X (i,j) m = λ m S + P (i,j) m (19)
A (i,j) m
is thus a N ×N -2 matrix obtained from removing the two columns corresponding to the sensors located at nodes (i, j).

X (i,j) m is vector X m without components i & j. Finally P (i,j) m = φm (T )X(T ) -φm (0)X(0) + A m . î X i + A m . ĵ X j (20)
where î & ĵ are the unit vectors having components îk = δ ki & ĵk = δ kj with δ ki the Kronecker symbol. Theorem 2 established in [START_REF] Caputo | Inverse source problem in a forced network[END_REF] states the unicity of the source : Theorem 2. For m ∈ N , N × N matrix A m defined in [START_REF] Pinto | Locating the source of diffusion in large-scale networks[END_REF], and ξ = {i, j}, if

• ξ is strategic • λ, ϕ m = 0 •
All N -2×N -2 sub-matrix obtained from extracting the two columns i et j and any two lines

s 1 & s 2 of A m denoted A (i,j) s 1 ,s 2 are invertible.
The data collected in ξ = {i, j} permit to uniquely identify S and λ.

Considering [START_REF] Diaconis | The cutoff phenomenon in finite Markov chains[END_REF], let us describe how obtaining S. First, it is important to mention that this is an over-determined problem having N equations for N -2 unknown. Source identification is obtained from using

-A (i,j) m X (i,j) m = P (i,j) m ( 21
)
This system is tested many times, by each time removing any couple of distinct lines s 1 and s 2 and solving (from invertibility condition of A

(ij) s 1 ,s 2 ) -A (i,j) s 1 ,s 2 X (i,j) m = P (i,j) s1 ,s 2 (22)
If two distinct coupled of lines (s 1 , s 2 ) and (s 1 , s 2 ) give the same X (i,j) m , then problem ( 19) is equivalent to [START_REF]The University of Florida Sparse Matrix Collection[END_REF], so that, it means that the source term is absent from [START_REF] Diaconis | The cutoff phenomenon in finite Markov chains[END_REF]. The repeated line in (s 1 , s 2 ), (s 1 , s 2 ) gives the location of the source s, i.e s is precisely the repeated line. The detailed algorithm of the source localization proposed in [START_REF] Caputo | Inverse source problem in a forced network[END_REF] is given in Algorithm 1.

Finally, we do not provide further details about the temporal form of λ(t), since the algorithm provided in [START_REF] Caputo | Inverse source problem in a forced network[END_REF] is strongly sensitive to added noise on the detector signal. It is interesting to mention the following limits of [START_REF] Caputo | Inverse source problem in a forced network[END_REF] :

• Strategic set is limited to two sensors, whilst it is possible to exhibit some graph without possible strategic set built from two sensor location only.

•

The regression matrix R can be rank deficient (the condition for regression is N s (T e -T 0 )/∆T > 2N ).

•

No indication are given for choosing T e the total recorded time and which mode µ m is best for projection.

•

The proposed method is very sensitive to added noise.

•

The algorithmic cost is O(N 4 ). [START_REF] Caputo | Inverse source problem in a forced network[END_REF] uses the working hypothesis that every sub-matrix

Extension to an arbitrary number of sensors

A (ij) s 1 ,s 2 are invertible for every couple of line s 1 , s 2 ∈ V sup- pressed from A (ij)
m . This hypothesis can be weakens so as to provide a broader generality to the source identification results when considering an arbitrary number of sensors located in ξ ⊂ V for which the corresponding hypothesis that every sub-matrix A ξ s 1 ,s 2 is invertible become a strong assumption. In fact some numerical test with N s = 3 sensors are sufficient to exhibit non-invertible sub-matrix A ξ s 1 ,s 2 . Hence there is a need for a technical softening of Theorem 2 hypothesis. In fact, the invertibility of sub-matrix A ξ s 1 ,s 2 is not strictly necessary. Let us define X m |ξ ⊥ the restriction of X m over the nodes where there is no sensors on V , and P s1 ,s 2 |ξ ⊥ , the similar restriction of P s 1 ,s 2 , One needs to find a solution for X m |ξ ⊥ for the generalization of ( 22) to an arbitrary number of sensors, i.e,

-A ξ s 1 ,s 2 X m |ξ ⊥ = P s1 ,s 2 |ξ ⊥ . ( 23
)
The over-constrained solution of ( 23) is unique when

Rank(A ξ s 1 ,s 2 ) = N -N s .
With this condition, one can then compute

X m |ξ ⊥ = -[(A ξ s 1 ,s 2 ) • A ξ s 1 ,s 2 ] -1 A ξ s 1 ,s 2 P s1 ,s 2 |ξ ⊥ . ( 24
)
Otherwise, if Rank(A ξ s 1 ,s 2 ) = N -N s several possible solutions for X m |ξ ⊥ could lead to disregard the true source s, resulting in non-uniqueness of the method. We now precisely state and prove this result.

Algorithm 1 Caputo et al. [START_REF] Caputo | Inverse source problem in a forced network[END_REF] algorithm for source identification Require: Single component Graph G(V, E) and its adja-

cency matrix A G 1: Chose two sensors location i, j ∈ V 2: if (i, j) is a strategic set then 3: Get X i (t) and X j (t) for t ∈ [0, T ] 4: Construct the Graph Laplacian matrix ∆ G from A G 5: Compute the set of eigenvectors v n ∈ V & eigenval- ues ω n ∈ Ω for n ∈ [1, N ] of ∆ G 6:
Compute the 2N-vector Y(T ) ≡ [y(T ), ẏ(T )], for T > T 0 , from inversion of the linear regression matrix [START_REF] Scardoni | Centralities Based Analysis of Complex Networks[END_REF] for R and (60) for D 7:

R : Y = [R • R] -1 R D Cf
Then, evaluate X(T )

8:
for each s ∈ V do 9:

for each q ∈ V \ {s } do 10: Evaluate A (i,j) s ,q & P (i,j) s ,q Cf (22) for A (i,j) s ,q for P (i,j) s ,q 11: Compute X (i,j) s ,q = -[A (i,j) s ,q • A (i,j) s ,q ] -1 A (i,j) s ,q P (i,j) s ,q 12:
end for 13: for each q, q ∈ V \ {s } do

14: Compute E (i,j) s ,q,q = X (i,j) s ,q -X (i,j) s ,q 2
15:

end for [START_REF] Peng | Smartphone malware and its propagation modeling: A survey[END_REF]:

Estimate

E (i,j) s = q,q ∈V \{s } E (i,j) s ,q,q 17:
end for 18:

for each s ∈ V do 19:

Find

s sol such that min s ∈V E (i,j) s = E (i,j) s sol 20:
end for

) ∈ S , Rank(A ξ s 1 ,s 2 ) = N -card(ξ)
, then a single observation M (λ, S) uniquely identifies the source node s.

Proof

Let (λ (1) , S (1) ) & (λ (2) , S (2) ) be two possible events associated with two sources candidates s 1 & s 2 ∈ S . These events leads to solutions X (1) & X (2) . Let us define Z = X (2) -X (1) . Using definition (5) of the observation operator, M (λ (1) , S (1)

) = M (λ (2) , S (2) ) ⇒ ∀k ∈ ξ, ∀t ∈ (0, T ), Z k (t) = 0.
(25) Using lemma 2.2 of [START_REF] Caputo | Inverse source problem in a forced network[END_REF] leads to Z(T ) = Ż(T ) = 0. Furthermore, from a unique observation D at ξ nodes, then the projection of mode decomposition Z m of Z on ξ is zero, i.e, ∀k ∈ ξ, Z k = 0, so that

A ξ m Z|ξ ⊥ = λ (2) m S (2) -λ (1)
m S (1) .

Removing lines

s 1 & s 2 from A ξ m leads to A ξ s 1 ,s 2 Z|ξ ⊥ = 0. (27)
Since

s 1 &s 2 ∈ S , Rank A ξ s 1 ,s 2 = N -N s , Rank (A ξ s 1 ,s 2 ) • A ξ s 1 ,s 2 = N -N s and one gets Z|ξ ⊥ = 0. ( 28
)
Hence, (λ (2) m , S (2) ) = (λ (1) m , S (1) ),

so that, in particular,

S (2) = S (1) , (30)
and uniqueness follows.

Extension to boundary condition sources (BCevents)

New definitions

Let us now consider the conditions needed for BC-events to be identified. We need to define the generalization of strategic set for BC-events, that we name BC-strategic set.

In this case, we need to enlarge the frame so as to consider all possible s-directed graphs G s ∈ G BC obtained from undirected graph G. This first necessitates the direct generalization of the previous notion of strategic set, although, insufficient in this case, that we thus name Componentdiscriminant set for BC-event, and define as

Definition 7. (Component-discriminant set) A set K 1 ⊂ V is component-discriminant if ∀G v ∈ G BC , ∀ω ∈ Ω, ∃k ∈ K 1 / v ω|Gv k = 0.
Furthermore, we need to define a vector-discriminant set Definition 8.

(vector-discriminant set) A set K 2 ⊂ V is vector-discriminant if ∀G v1 &G v2 ∈ G BC , ∀ω ∈ Ω v1 ∩ Ω v2 , v ω|Gv 1 = v ω|Gv 2 ⇒ ∃k = k ∈ K 2 / v ω|Gv 1 k v ω|Gv 1 k = v ω|Gv 2 k v ω|Gv 2 k
. so that a BC-strategic set can be defined Definition 9. (BC-strategic set) A set K ⊂ V is BCstrategic if it is both component-discriminant and vectordiscriminant.

BC-event uniqueness conditions

For BC-events at candidate node s , the corresponding Laplacian matrix ∆ G s associated with G s is distinct from ∆ G which is associated with the graph G. Let us denote ∆ G\s the N -1 × N -1 sub-matrix built from the N × N matrix ∆ G by removing its s th column and its s th line.

Let us define C ∆ G \s as being the s th column of ∆ G without its component s . Then, we define ∆ G s built from ∆ G\s as

∆ G s =         ∆ G\s C ∆ G \s 0 0         , (31)
where, on the last line, 0 denotes the N -1 component zero vector. This choice is made so as to permit the directedgraph forcing of BC-event from dynamic system

        Ẍ\s Ẍs         +         ∆ G\s C ∆ s 0 . . . 0 0                 X \s X s         =         0 . . . 0 λ(t)         , (32)
so that the s -neighbor nodes have no influence on the forced node, which, on the contrary, can propagate the forcing wave into the graph. This is why the event forcing λ(t) drives the source node acceleration term Ẍs in (32). Unicity is secured by initial conditions X(0) et Ẋ(0). Let us now consider the eigenvalues and eigenvectors of matrix ∆ G s . First, let us consider the R N unit-vector

1 N = 1 √ N [1, 1, . . . , 1]
. From the structure of ∆ G s given in (31), one can easily find that 1 N is a zero-eigenvector of ∆ G s since the unit-vector 1 N -1 in R N -1 is also the zeroeigenvector of ∆ G\s . So that (the first) zero-eigenvector of

∆ G s is v 1|G s = 1 √ N [1, 1, . . . , 1] .
Then, let us consider the sub-matrix ∆ G\s of ∆ G s having an orthogonal eigenvectors base (since it is symmetric)

denoted v n |G\s with n ∈ [1, N -1]. Let us then built a new base for ∆ G s denoted v n|G s = [v n +1|G\s , 0] with n = n + 1∀n ∈ [1, N -1]. v n|G s form an orthogonal base of R N -1 . They also are eigenvectors of ∆ Gs . Moreover v 1|G s is linearly dependant with v n|G s with n ∈ [2, N ], so that Im(v 2|G s , .., v N|G s ) spans R N -1 .
Hence, the (nonorthogonal) eigenvector base of ∆ G s is known.

Lemma 3. (Extension of theorem 1)

If (λ, S) is part of admissible BC-events set, i.e (λ, S) ∈ A s , X(t) solution of (3) given at nodes of BC-strategic set ξ, i.e X| ξ (t), permit to find a unique X(T ).

Proof

Let us consider two solutions s 1 (t) & X s 2 (t) associated with two source candidates s 1 & s 2 . From (10) the differences between these solutions, reads

If those solutions provide de same detected signal on ξ, ∀k ∈ ξ, we have

0 = y s 2 1 (T ) + (t -T ) ẏs 2 1 (T ) v 1|G s 2 k (34) -y s 1 1 (T ) + (t -T ) ẏs 1 1 (T ) v 1|G s 1 k (35) + N n=2 y s 2 n (T ) cos ω G s 2 n (t -T ) v n|G s 2 k (36) - N n=2 y s 1 n (T ) cos ω G s 1 n (t -T ) v n|G s 1 k (37) + N n=2 y s 2 n (T ) ω G s 2 n sin ω G s 2 n (t -T ) v n|G s 2 k (38) - N n=2 y s 1 n (T ) ω G s 1 n sin ω G s 1 n (t -T ) v n|G s 1 k (39)
For τ ∈ R + integrating (39) over [T, T + τ ], and dividing by τ 2 whilst taking the limit τ → ∞, one gets

ẏs 1 1 (T )v 1|G s 1 k = ẏs 2 1 (T )v 1|G s 2 k (40)
But yet, in every case, the eigenvector associated to the nul eigenvalue is

(1, 1, • • • , 1) 1/ √ N , so that v 1|G s 1 = v 1|G s 1 = (1, 1, • • • , 1) 1/ √ N , and thus ẏs 1 1 (T ) = ẏs 2 1 (T ).
Furthermore, integrating (39) in interval [T, T + τ ], whilst later on dividing by τ in the limit τ → ∞, one gets

y s 1 1 (T )v 1|G s 1 k = y s 2 1 (T )v 1|G s 2 k (41)
Thus we have, y

s 1 1 (T ) = y s 2 1 (T ). So that, from (39), one is leaft with ∀k ∈ ξ 0 = N n=2 y s 2 n (T ) cos ω G s 2 n (t -T ) v n|G s 2 k (42) - N n=2 y s 1 n (T ) cos ω G s 1 n (t -T ) v n|G s 1 k (43) + N n=2 y s 2 n (T ) ω G s 2 n sin ω G s 2 n (t -T ) v n|G s 2 k (44) - N n=2 y s 1 n (T ) ω G s 1 n sin ω G s 1 n (t -T ) v n|G s 1 k (45) So be it n ∈ [2, N ], let us project over ω G s 1 n mode from multiplying (45) by cos ω G s 1 n (t -T ) . Then, integrating over [T, T + τ ],
whilst, later-on dividing by τ and taking the limit τ → ∞, two distinct cases are possible

   ∀k ∈ ξ, y s 1 n (T )v n |G s 1 k = 0 if ω G s 1 n / ∈ Λ(G s 2 ) ∀k ∈ ξ, y s 1 n (T )v n |G s 1 k = y s 2 n (T )v n |G s 2 k if ω G s 1 n ∈ Λ(G 2 )
(46) The first case imply that y n (T )

s 1 = 0 because k ∈ ξ which is composant-discriminant, i.e ∃k such that v n |G s 1 k = 0. Furthermore, since ω G s 2 n / ∈ Λ(G s 1 )
, one also finds that y n (T )

G s 2 = 0.

The second case implies that ∃k = k ∈ ξ such that

v ω|G s 1 k v ω|G s 1 k = v ω|G s 2 k v ω|G s 2 k . But yet, if one supposes that v ω|G s 1 = v ω|G s 2 this implies ∃k = k ∈ ξ, being vector-discriminant, so that v ω|G s 1 k v ω|G s 1 k = v ω|G s 2 k v ω|G s 2 k
. Thus, reductio ad impossibile,

v ω|G s 1 = v ω|G s 2 .
We thus conclude that y

s 1 n (T ) = y s 2 n (T ). Thus, finaly since for all n ∈ [1, N ] one gets      y s 1 n (T ) = y s 2 n (T ) = 0, or y s 1 n (T ) = y s 2 n (T ) & v n|G s 1 = v n|G s 2 . ( 47
)
Unicity follows, i.e

X s 2 (T ) -X s 1 (T ) = N n=1 y n (T ) v n|G s 2 -v n|G s 1 = 0.
(48)

ALGORITHMS FOR GENERAL SOURCE IDENTIFI-CATION

Building general strategic sets

Here we describe the algorithm for finding a robust strategic set, both including the one defined in [START_REF] Caputo | Inverse source problem in a forced network[END_REF] for intern-type events, and BC-strategic sets (Cf Def 9) associated with BCevents. Algorithm 2 details a general method to built this strategic set.

Algorithm 2 Find Strategic set ξ for general sources and boundary conditions

Require: Single component Graph G(V, E) and its adja-

cency matrix A G Construct the Graph Laplacian matrix ∆ G from A G 2: for each v ∈ V do Construct the matrix ∆ Gv 4: Compute eigenvectors v n v ∈ V v & eigenvalues ω n v ∈ Ω v for n ∈ [1, N ] of ∆ Gv end for 6: for each v 1 , v 2 ∈ V do if Ω v1 ∩ Ω v2 = {∅} then 8: if (k, k ) ∈ ξ such that v n v 1 ,k v n v 1 ,k = v n v 2 ,k v n v 2 ,k then Find k = k ∈ V such that v n v 1 ,k v n v 1 ,k = v n v 2 ,k v n v 2 ,k 10: ξ = ξ ∪ {k, k } end if 12:
end if end for

General intern-source or BC-source identification

We now describe how the identification of general source of unknown type (intern or BC type) can be done. The identification method proceed from two main steps. A precomputational one associated with

•

Graph properties evaluation (Adjacency matrix, Laplacian eigenvectors)

• Regression computation to find modes amplitudes [Y(T ), Ẏ(T )] in [START_REF] Caputo | Oscillations of networks: The role of soft nodes[END_REF].

Then, a the source term identification step is performed where the prediction over the state X(t) from various source candidates is compared with the one observed on detectors locations. Algorithm 3 provides the details of the source identification implementation for unknown type source to be detected. This algorithm branch into the previous algorithm 1 of [START_REF] Caputo | Inverse source problem in a forced network[END_REF], in the special case of an identified intern type source. It also necessitates the computation of all possible Laplacian ∆ G s of every possible s -directed graphs G s .

EXTENSIONS AND IMPROVMENT FOR SOURCE IDENTIFICATION EFFICIENCY

Noise sensitivity

This section specifically consider the noise sensitivity of the source detection method.

Noise sensitivity of linear regression

As previously explained in §3.2, the detection method proceed from two main steps. The (pre-computational) regression step, and later on, the source term identification step. Let us first analyse the noise sensitivity of the regression step. We consider that the detector signal decomposes into Algorithm 3 General internal source or Boundary Condition (BC) identification Require: Single component Graph G(V, E) and its adjacency matrix A G Construct the Graph Laplacian matrix ∆ G from A G 2: Pre-sectorise V to get candidates set S ⊂ V Cf Alg. [START_REF] Giudicianni | Topological Taxonomy of Water Distribution Networks[END_REF] for pre-sectorisation for each candidate s ∈ S do return s sol , pertub_type_id end if Even if this computation is much more costly than the single graph Laplacian computation ∆ G needed for interntype source, it only represents a pre-treatment cost of the algorithm, which might also possibly be stored, so that it is not needed each time an identification is performed. Furthermore, s source candidates do not need to span the entire set V , but only a part of it, i.e S ⊂ V , a presectorisation issue later-on discussed in §4.3. the true signal D T plus some random noise vector η (associated with the discrete sampling of a continuous random process η(t)

4: Construct matrix ∆ G s Cf (4) for ∆ G s Compute eigenvectors v n v ∈ V v & eigenvalues ω n v ∈ Ω v for n ∈ [1, N ] of ∆ Gv 6: Estimate R s Compute Y s = [R s • R s ] -1 R s D with 2N- vector Y s ≡ [y s (T ), ẏs (T )] 8: Estimate R s • Y s -D 2 = r s end for 10: Find s sol such that min s ∈S r s = r s sol Compute Y = [R • R] -1 R D with 2N-vector Y ≡ [y(T ), ẏ(T )] 12: Estimate R • Y -D 2 = r if r <
D = D T + η (49)
Using the linearity of relation ( 13) leads to

Y = Y T + δY = [R • R] -1 R D T + [R • R] -1 R η (50) i.e δY = [R • R] -1 R η (51)
This extra-term δY is then responsible for an extraprediction δX, which, from the linearity of X with Y in [START_REF] Wang | An algorithm for locating propagation source in complex networks[END_REF] reads

δX(t) = δy 1 (T ) + (t -T )δ ẏ1 (T ) v 1 + (52) N n=2 δy n (T ) cos(ω n (t -T )) + δ ẏn (T ) ω n sin(ω n (t -T )) v n
We now examine how the extra-term δX(t) affects the second identification step of the algorithm.

Noise impact of source identification

The next step involves the projection of the state vector X(t) into the chosen mode m, denoted X m . Considering [START_REF] Diaconis | The cutoff phenomenon in finite Markov chains[END_REF] associated with the restriction of X m over non-detector nodes denoted

X m |ξ ⊥ , -A ξ m X m |ξ ⊥ = λ m S + P m |ξ ⊥ , (53)
One can realize that the right-hand-side term P m |ξ ⊥ also fulfills an extra-noise term that we now wish to evaluate. From using the generalization of ( 20) to an arbitrary number of sensors

P m |ξ ⊥ = φm (T )X(T ) -φm (0)X(0) + k∈ξ A m . k X k , ( 54
)
one is able to compute the extra-noise contribution of 24) to an extra-prediction given by :

P m |ξ ⊥ δP m |ξ ⊥ = φm (T )δX(T ) -φm (0)δX(0) + k∈ξ A m . k δX k , (
δX m |ξ ⊥ = -[(A ξ s 1 ,s 2 ) • A ξ s 1 ,s 2 ] -1 A ξ s 1 ,s 2 δP s1 ,s 2 |ξ ⊥ (56)
This extra-prediction is responsible for the noise sensitivity of the detection algorithm. ( 56) and (55) combined with (52) and finally (51) provide the linear relationship between this extra-prediction δX m |ξ ⊥ and the noise vector η, i. e. the sensitivity matrix of the state prediction over noise. Since, this sensitivity matrix involves the X k , we now study how its sensitivity can be attenuated when increasing the timeprojection support T , or conversely increasing the recording frequency.

Sampling frequency and detectors number effect

on X k and δX k (T )

Sampling frequency effect

First, let us note that the noise is weaken by the scalar product for

δX k = T 0 φ m (t)δX k (t)dt, (57)
with, as mentioned just above, δX k (t) linearly related to vector η(t). For a centered, normalized, decorrelated random noise process η k (t), i.e

Cov(η k (t)η k (t )) = E(η k (t)η k (t )) = σ 2 N δ(t -t )
, (58) the central limit theorem provides an asymptotic limit for the behaviour of δX k as T → ∞ such as E(δX k 2 ) ∼ 1/T .

Similarly, for fixed T when increasing the recording frequency f , one expects E(δX k 2 ) ∼ 1/f . We numerically test this prediction using numerical computations on a card(V ) = N = 10 nodes random graph, having maximum vertex degree ∆(G) = 3, minimum vertex degree δ(G) = 2, with two, randomly located, detectors in figure 3a (for k = 1). It is interesting to observed how well the expected asymptotic rapidly settles for finite sampling rate as small as f = 20 Hz for a finite time-support T = 100s and a final time detector recording T e = 120s, and a fixed noise r.m.s amplitude σ N = 10 -2 . We furthermore analyze the sensitivity of δX k (T ) on the sampling frequency which is similarly ruled by a similar central-limit asymptotic behavior as shown in figure 3b.

Noise sensitivity versus detector number

For a fixed σ N = 10 -2 the increase of detector number leads to an increase of sampling data in the regression steps. This leads to a decrease of the apparent noise of δX(T ). Nevertheless, this decrease falls short after the detector number reaches a (small) fraction of the total node number N , where after the benefits of increasing their number almost saturates. This effect is illustrated in tests on randomly generated graphs indicates that, after the detector number has reached a small fraction of N , there is a a very modest interest for increasing further their number. More precisely we numerically found a sharp transition observed in Figure 4a nearby N s N/12 for which the decrease of the induced noise becomes very weakly sensitive to a further increase in detector number. We believe this can be understood as a statistical effect associated with a partial correlation of observed signals when detector spatial density increases. Since recorded signal are used in the regression process if observations are correlated, their noise becomes weakly sensitive to time averaging. This effect is reminiscent of the cut-off phenomenon arising in Markov chains shuffling [START_REF] Diaconis | The cutoff phenomenon in finite Markov chains[END_REF] and might be directly connected with it. If so, the detector cut-off should scale logarithmically with N . In any case, we believe this cut-off in detector de-noising efficiency to happen in any graph, for any detector distribution. This is an important finding for the design of optimal strategies for detector number choice in a given network. This now leads us to explore the effect of noise r.m.s on the algorithmic efficiency.

Designing detector number versus noise amplitude

At finite sampling rate f , for finite time-support T = 100s and final time detector recording T e = 120s, we analyze the effect of the noise r.m.s σ N on the success rate of the detection algorithm. Indeed, when the noise r.m.s increases, the source detection is not always possible since the signal is degraded. Furthermore, in real applications the noise r.m.s is rarely smaller than 10 -2 so that it is relevant to analyze how the method performs when increasing noise. Figure 5 provides this analysis for a sampling rate lower than 50Hz (chosen as the relevant maximal sampling rate) for three detector relative number N s /N = 10%, 15% and 20%. The chosen initial condition is λ(t) = 2∂ 2 t tanh((t -40)/10) where 40s is the chosen trigered time for the event, having a typical time variation of 0.1s. This figure showcase an abacus of practical significance, since for a given noise amplitude, providing the number of sensors (and the other sampling and averaging parameters T & T e ) one can design the sampling frequency for which 95% of success identification is reached. It is interesting to note that the progress obtained from having N s = 15 to N s = 20 is rather poor. The propagating time is then used backward to evaluate a distance on the graph We believe this observation related to the effect previously illustrated in figure 4, i.e a sharp transition in sampling efficiency when the detector relative number is larger than 8%. Figure 5 confirms that there is a cutt-off density efficiency for detector's noise reduction improvement.

Presectorisation

General principle

The general idea of pre-sectorisation is to transpose triangulation methods used in Euclidian space into graphs. This suppose that the propagating wave inside the network indeed spend physical time along each edge traveling from one node to reach its neighbor. This propagating time is used to set a distance to each edge of the graph enforcing a metric onto the graph. as proposed in [START_REF] Shen | Locating the source of diffusion in complex networks by time-reversal backward spreading[END_REF]. From each detected signal at any given detector, and provided the time origin of the unknown event, one is able to evaluate a wave-front of possible source positions. The intersection of these wave-fronts obtained from several detectors able the source pre-sectorisation (Cf Figure 6a). Nevertheless since the time origin t 0 of each event is unknown, one needs to set t 0 as an unknown parameter of the pre-sectorisation method. This result in possible difficulties exemplified in Figure 6b,c,d. Since the traveling time permits to define a metric into the graph, one can built the convex-hull of the detectors set inside the graph. Considering this convex-hull, two possibilities arise. Either, the event appears inside this convex-hull (Cf Figure 6a,b,c), or outside it (Cf Figure 6d,e). In the first case, the retrograde wave-fronts intersection, emergent from each detectors, are spreading along a limited area, resulting in a moderate number of source candidates. In the second case, if the convex hull only represent a moderate fraction of the graph, the source being outside it, almost any node outside the convex hull can provide a possible source, resulting in a large number of candidates. Hence, for an efficient and robust pre-sectorisation, it is important that the convex-hull of the detector set covers a large portion of the entire graph. This is a supplementary and practical constraint on the detector set, additional to the mathematical requirement that is should also be BCstrategic. This is why considering a detector number larger than two is a relevant practical issue.

We now describe in more details the chosen presectorisation method.

Method and algorithm

As previously mentioned, the event's time origin t 0 is unknown when the first arrival signal is detected. Let us define t 0g for this time origin guess. Let us call S as the set of all possible candidates, and S (t 0g ) the set of all possible candidates obtained from pre-sectorisation whilst supposing the time origin is t 0g . Then, one has

S = tmax t0g=0 S (t 0g ). (59)
To find S (t 0g ), a retrograde wave-fronts propagation is computed from every detector in ξ. Let us define d k (t 0g ) the distance between the source and each detector k ∈ ξ and the distance matrix D between any node v and detector k.

If the graph is embedded in space (such as a pipe network for which the vertex to vertex distance along the graph is a 'natural' feature of the graph), this distance is the "Euclidian/topological" distance: i.e the sum of the Euclidian distances from every successive couple of adjacent nodes along the shortest path into the graph between node n and detector k. If the graph is not embedded in space it can be defined as the traveling wave distance define over the (ordered) shortest path S p joining every successive vertex pairs v , v along itself and associated with a traveling velocity c v ,v at each edge (in most of the formulation provided in this paper we have supposed an homogeneous traveling wave velocity c, i.e c v ,v = c, as well as an homogeneous edges, i.e a non-weighted graphs, but the generalization to weighted graph is direct, this is why, we mention it here) and time traveling ∆T v ,v along each edges from vertex v to v

D ≡ D kv = v ,v ∈Sp(k→v) c v ,v ∆T v ,v . (60)
Then one can define S (t 0g )

S (t 0g ) = k∈ξ {v ∈ V \ D kv = d k (t 0g )} (61)
Albeit theoretically sound, this expression does not include some tolerance on the signal detection due to signal noise or detector sensitivity. Including an extra parameter ε, one can define a more physically relevant version of (61)

S (t 0g ) = v∈ξ {v ∈ V \ |D kv -d k (t 0g )| < ε)} (62)
So that, (59) with (62) provide a robust and fairly general (independent of specific graph features) formulation so as to locate the source into a pre-sectorised region. Other more dedicated, specific and more effective methods could possibly be set-up to achieve a similar pre-sectorisation. Nevertheless, the one presented here, already provide a drastic improvement over a systematic exploration algorithm, since it is easy to show that it permits a O(N ) gain in complexity on the source detection algorithm over the systematic exploration proposed in [START_REF] Caputo | Inverse source problem in a forced network[END_REF]. First let us summary the complexity of a systematic exploration when the presectorisation set is the entire graph, in algorithm 3. In this case in step 2, the subset of source candidates spans the entire vertex set, i.e S = V so that the loop over each candidates is O(N ). The subsequent steps of algorithm 3, inside this loop, require

Step 4. Construction of matrix ∆ G s : O(N 2 ),

Step 5. Eigenvectors v n v computation :

O(N 3 ), Step 6. Building matrix R s : O(N.N s ), Step 7. Compute Y s = [R s • R s ] -1 R s D : O(N 3 ), Step 8. Estimate R s • Y s -D 2 = r s : O(N.N s ),
Step 10. Find s sol such that min s ∈S r s :

O(N ), Step 11. Compute Y = [R • R] -1 R D : O(N 3 ), Step 12. Estimate R • Y -D 2 = r : O(N.N s ).
Thus the overall complexity from Step 4 to Step 12 is O(N 3 ). Including Step 4 to Step 12 in a O(N ) loop results in a O(N 4 ) cost associated with a systematic source candidates exploration, as found in [START_REF] Caputo | Inverse source problem in a forced network[END_REF]. Furthermore, since the pre-sectorisation restrict the source candidate to a finite number of nodes the size of which being independent of N , the external loop complexity over pre-sectorised candidates is O(1) leading to an overall O(N 3 ) cost of the pre-sectorisation algorithm. To be exhaustive, one should add the complexity of the pre-sectorisation computation being a pre-requisite of algorithm 3. The algorithm for presectorisation is hereby given in algorithm 4. resulting in an overall O(N s .N ) cost of this pre-computation step. Hence, as expected, the pre-sectorisation has a benign consequence on computational complexity.

Numerical results and illustrations on graph models

This section examines the algorithmic cost and efficiency on various graph models. The effect of pre-sectorisation on the computational time for source detection is illustrated in figure 7 for poorly connected random graphs. This Fig. 7: Log-log plot of computational time versus nodes number N . Continuous curves display the ∼ N 4 (purple) and ∼ N 3 (red) behaviors. Square symbols represent the average results over 10 random sources distributed within a random graph, with maximum vertex degree ∆(G) = 3, minimum vertex degree δ(G) = 2, with σ N = 5.10 -4 having N nodes for the pre-sectorisation method (error bar here indicate the min and max time over ten simulations). Bullets symbols represent the result of the systematic exploration algorithm of [START_REF] Caputo | Inverse source problem in a forced network[END_REF]. figure compares a systematic exploration of the network for finding the source used in [START_REF] Caputo | Inverse source problem in a forced network[END_REF] to the pre-sectorisation method proposed in the previous section. A significant number of randomly chosen source location are explored so as to provide an average numerical cost. Figure 7 convincingly illustrates that a O(N 3 ) computation cost can be achieved for large N (here the maximal value of N was 2000), as opposed to the systematic exploration cost being O(N 4 ). More precisely, we found that, for the chosen weakly connected graphs (having in mind pipe networks applications) the observed efficiency is somehow better than these theoretical predictions. Furthermore the success rate of the method is analyzed when varying the sensor density choosing an homogeneous random distribution of sensors within random graphs tested in [START_REF] Shen | Locating the source of diffusion in complex networks by time-reversal backward spreading[END_REF] : Erd ös-Rényi (ER) and a Watts-Strogatz (WS). It is found that the success rate of location events display a very similar behavior in both cases. More precisely, when varying the average node degree k (average edge per node) (using a probability p = 0.1 for WS graphs), Figure 8a indicates more than 80% of success for sensor density being larger than 15%. For sensor density smaller that 15% we found a rather sharp degradation of the success rate. However this conclusion is tempered by the observation of the "almost-found" (one degree distant along the graph) sources identification rate. The later is more robust to sensor density, being larger than 90% for sensor density equal to 10%. A similar sharp transition at lower sensor density is also observed for "almost-found" success rate. A similar conclusion can be reached for the Erd ös-Rényi graph family in Figure 8b, for which sensor density equal to 10% permits 80% of success in most cases for "almost-found" sources, whereas a 15% sensor density is needed to reach this achievement for exact findings. These performances of our method are close to the ones reported in [START_REF] Shen | Locating the source of diffusion in complex networks by time-reversal backward spreading[END_REF] for the same test cases.

APPLICATION TO REAL GRAPHS

The ability of the presented method to identify randomly chosen events within real graphs is illustrated in this section. Concerning real graphs, power networks as well as river basin networks (as in [START_REF] Pinto | Locating the source of diffusion in large-scale networks[END_REF]) taken from open data bases [START_REF] Davis | The University of Florida Sparse Matrix Collection[END_REF], have been specifically chosen so as to analyze two complementary class of graphs : a first class with a high density of circuits (power-network graphs) and a second one almost free of circuits having three-like structure (river basin networks). For every case, a sensor r.m.s noise equal to σ N = 1.5 10 -2 is chosen in order to compare the method's performances in distinct networks having similar noise level. These two classes of graphs also illustrate two possible classes of applications of the method. The ability to detect BC-events is evaluated here. Figures 9a,b,c,d thus illustrate three power-network distribution and one urban connectivity network (more details on graphs are given table 1) whereby the BC-event detection method with presectorisation has been applied. In each case, 0.15 × N sensor node number (a sensor density of 15%) have been used. The node number of the chosen real graph varies between N = 47 to N = 443 (as compared with N = 287 examined in [START_REF] Pinto | Locating the source of diffusion in large-scale networks[END_REF]).

In almost every cases, the source is exactly recovered, (100% success in Figure 9b,c,d,e,f). Moreover, location failure, if any, are always found very close to the exact source, as found in Figure 9a where only one degree distance is observed between the source and its false prediction. The dependence of success-rate with sensor density is more accurately investigated in Fig. 10, averaging each rate over N/2 randomly seeded sources. Fig. 10 shows that 10% sensor density already permits more than 90% success rate in every case.

CONCLUSION

This paper presents various extensions of event identification for non-dissipative wave propagation in networks described by graphs. Having defined intern events and BCevents we found that in both cases, for an arbitrary number of detectors 1) Source identification of intern-type events (resp. BCevent) are possible (and unique) if the detector set is strategic (resp. BC-strategic). This result extents the analysis of [START_REF] Caputo | Inverse source problem in a forced network[END_REF] dedicated to intern-type events only recorded by two sensors only. 2) We found that the noise sensitivity of the source identification detection can be improved in various ways. Increasing the time-averaging of the projection support, the time-recording of the detectors, or the frequency acquisition are all beneficial, with central-limit asymptotic behavior extending the analysis of [START_REF] Caputo | Inverse source problem in a forced network[END_REF]. We also estimate the sensitivity matrix response of noisy signals over an arbitrary set of detectors, an analysis not provided in [START_REF] Caputo | Inverse source problem in a forced network[END_REF].

Another original result obtained in this study has shown that increasing the detector number improves the noise robustness of the detection, up to a cut-off detector density being a small fraction of N (possibly scaling as ln(N )). 3) We analyze a pre-sectorisation strategy so as to improve the computational cost of source identification. We have stressed the interest of a uniform coverage of detectors over the network. 4) Finally, using pre-sectorisation we found that it is possible to improve the source detection algorithm within a card(V ) = N network to an O(N 3 ) computation cost, as previously obtained in the context of diffusive processes in graphs [START_REF] Shen | Locating the source of diffusion in complex networks by time-reversal backward spreading[END_REF], [START_REF] Pinto | Locating the source of diffusion in large-scale networks[END_REF]. The proposed method has been tested in real-world networks for which the applicability of the method on noisy data has been illustrated. Finally it is important to mention that, since the presented method lies upon a time projection upon a discrete base decomposition it necessitates the choice of a particular projection mode. The choice of this mode can affect the efficiency of the method.

Fig. 1 :

 1 Fig. 1: Illustration of graphs resulting from intern type or boundary condition type sources. (a) Non-oriented unchanged graph G with inner source node #3 (red) as in [1]. (b) Oriented graph resulting from a boundary condition source (BC-event) node #3 (red). Notation G 3 is introduced to describe this locally oriented graph at node #3 distinct from graph G.

Definition 4 (Fig. 2 :

 42 Fig. 2: Illustration of non-detectability of some sources resulting from the graph symmetries. (a) For sensor in node #3 (in black), a source in node #2 (b) or #4 (c) (in red) produce the same signal. Hence node #3 only does not provide a strategic set, i.e ξ = {3} is not strategic.

  55) where obviously, δX m ≡ δX, ϕ m with δX k the kth component of δX m . Hence during the identification procedure associated with the inversion of the linear problem (53), this extra term in δP leads, from (

Fig. 3 : 1 )

 31 Fig. 3: Root mean square of the variance estimator of X 1 , i.e E(δX 2 1 ) and δX(T ), i.e E(δX 2 1 (T )) versus the sampling frequency f for fixed σ N = 10 -2 , T 0 = 80s, T = 100s, T e = 120s, and N = 10. Dots are numerical results. Continuous line are the expected ∼ 1/ √ f asymptotic behavior.

Figure 4 Fig. 4 :

 44 Fig. 4: Impact of sensor number N s on the noise of X(T ), i.e the square root of the variance of δX(T ) for fixed σ N = 10 -2 , T 0 = 80s, T = 100s, T e = 120s, f = 100Hz and N = 100. (a) shows the sharp transition nearby N s = 8 in semi-log representation. (b) Same data as in (a) but in normal plot within N s range between 10 to 95.

Fig. 5 :Fig. 6 :

 56 Fig. 5: Frequency sampling f providing a 95% success for source identification versus r.m.s noise σ N for randomly generated graphs having maximum vertex degree ∆(G) = 3, minimum vertex degree δ(G) = 2, N = 100 nodes and randomly chosen sources.Dotted lines are quadratic fit of the results to help delineate the upper region for which the f -σ SN values offer a 95% confidence for source detection.

Algorithm 4

 4 Pre-sectorisation to obtain candidate node set S Require: Single component Graph G(V, E) and its adjacency matrix A G Require: Parameter N s -vector trig_time, BC-strategicN S -vector set ξ, tol Evaluate N s × N matrix D = D i,j for i ∈ ξ, j ∈ V Cf (60) for definition of D Initialization S = NULL,S t0 = NULL 3: for t 0 ∈ [min trig_time-max D/c, min trig_time] by ∆t 0 do S t0 = ∩ k∈ξ {v ∈ V /D kv -c(trig_time[k] -t 0 ) < tol } S = S ∪ S t0 6: end forThe complexity of algorithm 4 is as follows :Step 2. Evaluation of N s × N matrix D : O(N s .N )Step 3. Back propagation loop over sensors set : O(N s .N ),

Fig. 8 :

 8 Fig.8: Success rate of source detection with sensors having r.m.s noise σ N = 1.5 10 -2 versus sensor density within model networks used in[START_REF] Shen | Locating the source of diffusion in complex networks by time-reversal backward spreading[END_REF]. Continuous curves provide the success for exact source location. Dotted curves the success rate for one degree distance of source location. N = 100 nodes have been chosen as in[START_REF] Shen | Locating the source of diffusion in complex networks by time-reversal backward spreading[END_REF], and N = 100 test events have been considered. (a) Watts-Strogatz (WS) graphs with average node degree k (using a probability p = 0.1 for supplementary random connections). (b) Same conventions as (a) for Erd ös-Rényi (ER) graphs.

Fig. 9 :

 9 Fig. 9: Illustration of real graphs distribution networks from public data bases. Each exactly detected source is depicted in red whereas wrongly detected ones are represented in purple with a link to the exact (red) source -Cf (a)-. Sensor positions are randomly chosen with a sensor node density equal to 15% and r.m.s noise level σ N = 1.5 10 -2 . (a,b,c) are three power distribution networks. (d) is a urban connectivity network. Total node number N and graph origin are indicated in each case. (a) [21]N = 118 (b) [22] N = 274 (c) [23] N = 443 (d) [24] N = 209 (e) Garonne river basin (South of France) [25] mapped with Geo-spatial coordinates [26] N = 293 (f) same as (e), Garonne river, in Toulouse town. [27] N = 47.

Fig. 10 :

 10 Fig. 10: Detection success rate versus sensor density within real networks using the same conventions as Fig. 8 tested over N/2 events. N = 188 corresponds to power network of Fig. 9a, N = 209 to urban network of Fig. 9d, N = 293 to river basin network Fig. 9e.

  r s sol then

	14:	pertub_type_id=intern source type
		Apply steps 8 to 21 of Alg. 1 for s ∈ S
	16: else
		pertub_type_id =BC
	18:

APPENDIX

Franck Plourabou é is Senior Scientist at CNRS (DR1). He received his degree from Ecole Normale Sup érieure in Paris in 1993, and his PhD in Physics at Paris 7 University in 1996. Since then, he works at the Fluid mechanics Institute of Toulouse (IMFT UMR 5502), being specialist of heterogeneous media, micro-hydrodynamics, multi-physics modeling and complex networks arising in biological micro-circulation and hydraulic systems. Matrix derived from the graph Laplacian matrix acting on V -space so that card(Im(Am)

Un-directed Graph defined by Vertex and Edges sets

Compact support source perturbation temporal variation Strategic nodes set where signal is detected,

Restriction of X(t) at sensor's location ξ ⊂ V , dim(X|ξ) = Ns X|ξ ⊥ Restriction of X(t) out of sensor's location ξ ⊥ = V \ ξ X (i,j) Restriction of X(t) out of sensor's location (i, j) y(t)

Mode projection vector of


domain: info shs

text:
Luke Dzwonczyk

Léo Chédin

Alejandro Saldarriaga-Fuertes

Max Sherr

Hélène-Camille Crayencour

Carmine-Emanuele Cella

Source Separation Methods for Computer-assisted Orchestration

In this paper, we study the possibility of adding source separation as a preprocessing step to the computer-assisted orchestration process. We first discuss the motivation of this addition and its potential to increase the quality of orchestrations of multi-layered sounds. Second, we select several state-of-the-art models for both music source separation (separation of instruments) and universal sound separation (separation of arbitrary sounds), and compare their effectiveness for the task of orchestration. We assess which methods best suit the needs of orchestration by applying them on our own target sounds, orchestrating the separated outputs, and finally comparing them to the orchestration of the same target without separation. Our experiments show that source separation improves the quality of orchestrations, and the more accurate the separation, the better the resulting orchestration. Finally, we compare unsupervised methods to supervised methods for separation, and comment on the effect of training data selection on performance of supervised methods.

Introduction

Computer-assisted composition is a field that focuses on the creation of computational tools to aid in the musical composition process [START_REF] Fernandez | AI methods in algorithmic composition: A comprehensive survey[END_REF][START_REF] Ariza | Navigating the landscape of computer aided algorithmic composition systems: a definition, seven descriptors, and a lexicon of systems and research[END_REF]. Within this field lies target-based computer-assisted musical orchestration, a method which uses orchestral instruments to recreate the timbre and temporal evolution of a target sound. It is the process of creating an orchestral score that best matches an arbitrary target sound given a similarity metric and constraints [START_REF] Maresz | On computer-assisted orchestration[END_REF]. In other words, finding which instruments playing which notes sounds the most similar to the given sound. A goal of computer-assisted orchestration is to help composers by accelerating certain creative processes. In particular, target-based computer-assisted orchestration helps composers explore new timbral possibilities by combining instrumental samples in such a way as to mimic the timbre of a given target sound.

Target-based computer-assisted orchestration is a complex problem because it attempts to recreate the timbre of different sounds, and timbre is a difficult parameter to model. The task becomes more complicated when you consider multi-layered sounds that are a combination of multiple sources.

3rd Conference on AI Music Creativity (AIMC 2022) When two or more sources overlap in time, the orchestration algorithm cannot distinguish between them. As such, the orchestration is rigid: the whole orchestra suddenly changes configuration when a new source dominates, without any blending with the previous dominating source.

Our hypothesis is that if source separation is applied as a pre-processing step, then the separated layers of a sound can be independently orchestrated and then recombined, resulting in an improved orchestration. This logic follows from the orchestration principle of dovetailing (Adler, 2016, p. 467-473). Dovetailing is a well known practice aimed at creating a blend of the instruments used during orchestration in order to transition between timbres. Technically, dovetailing is the swapping and overlapping of musical lines between different instruments: musical lines are scattered across multiple instruments, connected with overlapping pitches Orchidea is a framework and set of tools that is currently considered the state-of-the-art for computerassisted orchestration. It embeds the target in a feature space and uses a jointly-optimized heuristic and constraint solver to find a combination of samples that best match the target [START_REF] Cella | Open-source modular toolbox for computer-aided orchestration[END_REF][START_REF] Cella | Orchidea: a comprehensive framework for target-based assisted orchestration[END_REF]. All of our orchestrations are done using Orchidea's command line tools for batch processing, which is available for download at http://www.orch-idea.org.

Orchidea is capable of performing both static and dynamic orchestration. In static orchestration, temporal evolution is ignored and the target is represented as a single vector of timbral features. The output is a single onset of notes that does not change over time. In contrast, dynamic orchestration considers a time series of features, which allows the orchestration to evolve over time and contain multiple onsets of notes. Orchidea is able to perform dynamic orchestration by cutting the target in time into multiple segments, and then orchestrating each segment. In this paper, we perform only static orchestrations in order to better observe the effects of source separation. The reason for this is that we wish to disentangle the two difficult problems of source separation and time segmentation. By only performing static orchestrations, we reduce the complexity of the problem by removing time segmentation and can better observe the effects of source separation.

The paper is organized as follows: in section 2 we introduce sound source separation and list other musical tasks in which source separation has had success. Section 3 will discuss our methods and walk through the testing procedures we employed in our experiments. Sections 4 and 5 will describe the findings of our experiments, an analysis of the results, and offer possibilities for future work.

The code for this paper can be found at: https://github.com/dzluke/AIMC2022, and you can listen to a selection of targets and orchestrations here: https://dzluke.github.io/AIMC2022/.

Background

Sound source separation is the process by which a single audio file is separated into multiple sound sources. A perfect separation is able to dissect a sound exactly into its constituent parts, in which each part is an independent sonic event. Music source separation is a specific application of source separation in which the input audio is music that is comprised of a subset of instruments or sounds. For example, one music source separation method could attempt to split a song into three parts: voice, bass, and drums. Another aimed at orchestral music may try to separate the input into families of instruments: woodwinds, brass, strings, and percussion. In contrast, universal sound separation does not operate under the assumption that the input is of a musical nature, and can be applied to arbitrary sounds. In this case, the number of sources expected is specified, and the method will attempt to divide the input into the given number of sources.

Using source separation techniques as a pre-processing step has been shown to improve results in several music processing tasks. Source separation can be used to remove part of the spectrum in order to enhance the clarity of the features of interest for the task at hand. For instance, in the task of automatic chord estimation from audio, some authors have chosen to remove part of the spectrum related to percussive sounds to improve chord accuracy [START_REF] Reed | Minimum classification error training to improve isolated chord recognition[END_REF]. Similarly some authors have investigated how beat tracking can be improved by using source separation as a pre-processing step for difficult songs with highly expressive vocals [START_REF] Zapata | Improving beat tracking in the presence of highly predominant vocals using source separation techniques: Preliminary study[END_REF]. Source separation can also be applied to separate different components of the signal and work on them separately to reduce the complexity of the task. For instance, in the context of tempo detection in [START_REF] Chordia | Using source separation to improve tempo detection[END_REF], the authors decompose the signal into sources to reduce rhythmic complexity, under the idea that some layers may be more rhythmically regular than the overall mix. Another example can be found in [START_REF] Gómez | Jazz solo instrument classification with convolutional neural networks, source separation, and transfer learning[END_REF] where harmonic/percussive and solo/accompaniment source separation techniques are investigated as a pre-processing step to improve predominant instrument recognition in jazz music.

Consider a motivating example: a target that has a continuous drone sound throughout and a melody playing above this drone. Without source separation, Orchidea would detect each note of the melody as a new onset and cut the target in time at each new note. However, this segmentation in time would also affect the drone, cutting it in time and forcing it to have a new onset each time the melody changes notes. The resulting orchestration would have multiple onsets for the drone, even though only one onset is needed at the beginning of the drone. When source separation is applied to this target, the drone and melody could be orchestrated separately. Therefore, the drone could be orchestrated by a single onset that lasts for the duration of the drone. At the same time, the melody could have as many new onsets as needed without affecting the drone. The orchestration would be greatly improved by removing the unnecessary onsets in the orchestration of the drone. We compare the effectiveness of different source separation methods for the task of orchestration by applying the separation methods to various targets, orchestrating the separated output of the method, and finally comparing this to the orchestration of the target if no separation was performed. This process is described in Fig. 1. Consider a target sound that is a combination of a low droning sound and high-pitched whistle. With separation, these two sounds would be disentangled into two separate sub-targets and each would be individually orchestrated. The separated orchestrations would then be combined to create the solution (Fig. 1,bottom). This would then be compared to the orchestration if no separation had been performed (Fig. 1,top).

Experimental Methodology

We test the effectiveness of five different supervised source separation methods on 300 custom made targets that are combinations of sounds from different databases. During testing, a target is input to a source separation method, which outputs four sub-targets. The accuracy of the separation is measured using the cosine distance. Then each sub-target is independently orchestrated with a randomly assigned subset of the full orchestra. These orchestrations are then combined to play simultaneously, creating a final orchestrated solution. Finally the distance between the target and solution is calculated, giving us a metric to compare the various separation methods. See Fig. 2 for a diagram of this process.

For a given target, the same split of the orchestra is used to orchestrate the separations. We randomly create these subsets of the orchestra in order to avoid any biases that could come from a hand-picking the instruments used to orchestrate each sub-target.

Data

We created 300 targets where each target was a combination of four source sounds. The sources come from the NIGENS [START_REF] Trowitzsch | The NIGENS general sound events database[END_REF] and BBC (BBC, 2021) databases, and freesound.org [START_REF] Font | Freesound technical demo[END_REF]. We selected a total of 90 samples from these databases, choosing sounds that fit the following criteria:

1. Static sounds that do not change harmonically over time 2. Sounds in which there is at least some pitched content and not only noise

The sounds chosen include alarms, bells, engine noises, soundscapes, synthesizer chords, and sound effects. Each target that we used for testing was a combination of four randomly chosen source sounds from the group of 90 sounds. When a target is created, the longest of the four sources is chosen to start playing at the beginning of the target. The other three sources are then randomly assigned different times to begin playing such that they will start between the beginning of the target and half-way through the longest source sound. This is done to ensure that the sources overlap and to minimize the amount of time in which there is only a single source sounding.

Separation methods

We consider source separation models trained to tackle two distinct separation problems: music source separation and universal sound separation. The music source separation models we test are Open-Unmix, Demucs, and Conv-Tasnet. The universal source separation models we test are non-negative matrix factorization (NMF) and TDCN++. Open-Unmix, Demucs, Conv-Tasnet, and TDCN++ are supervised models, and NMF is an unsupervised method.

It is important here to make the distinction that we did not train the supervised models ourselves. Therefore, we cannot guarantee that the data the models were trained on is optimal for solving our problem. However, one of the goals of this paper is to explore whether existing, pre-trained source separation models can separate non-musical sounds in a way that improves orchestration. We choose to include music source separation methods because of the diversity of different models available and the results these models achieve at their task.

Music source separation

The music source separation models we use are Open-Unmix [START_REF] Stöter | Open-Unmix -A Reference Implementation for Music Source Separation[END_REF], Demucs [START_REF] Défossez | Music source separation in the waveform domain[END_REF], and Conv-Tasnet [START_REF] Luo | Tasnet: Surpassing ideal time-frequency masking for speech separation[END_REF]. For these architectures, we use models that are pre-trained on the MUSDB18 dataset [START_REF] Rafii | The MUSDB18 corpus for music separation[END_REF]. The models output 4 stereo tracks that correspond to instrumental categories defined in the SiSEC 2018 Mus evaluation campaign: vocals, drums, bass and other [START_REF] Stöter | The 2018 signal separation evaluation campaign[END_REF].

Open-Unmix [START_REF] Stöter | Open-Unmix -A Reference Implementation for Music Source Separation[END_REF] is an open source implementation of a deep music source separation model, based on [START_REF] Uhlich | Improving music source separation based on deep neural networks through data augmentation and network blending[END_REF]. It is actually composed of multiple models that are trained for each instrumental target. Each of these models is trained on the specific target using the same architecture, based on a three layer bidirectional LSTM network. Open-Unmix operates on the Short Time Fourier Transform (STFT) of the input mixture. It predicts the target by multiplying the output of the LSTM network with the magnitude spectrogram of the input, essentially applying a mask, and uses Wiener filtering as a last processing step.

The Demucs [START_REF] Défossez | Music source separation in the waveform domain[END_REF] architecture is a deep learning model for musical source separation. This models in the waveform domain, taking the stereo mixture as input, and outputting 4 stereo waveforms corresponding to the four categories of the MUSDB18 dataset. The Demucs models consists of a 6-layer convolutional encoder, two bidirectional-LSTM layers, and a 6-layer decoder made of transposed convolution. Skip U-net connections [START_REF] Ronneberger | U-net: Convolutional networks for biomedical image segmentation[END_REF] link encoder layers with their corresponding decoder layers. This model obtained state-of-the-art results on the MusDB set, surpassing Open-Unmix. We use the second version of Demucs which uses only the waveform, the third version being a hybrid model using both the waveform and a spectrogram.

We also used the Conv-Tasnet (Luo and Mesgarani, 2018) implementation of [START_REF] Défossez | Music source separation in the waveform domain[END_REF], who adapted the architecture, originally designed for speech separation at 8 kHz, for the task of music source separation. Conv-Tasnet is a deep learning model that is also based on three processing stages: encoder, separation and decoder. The encoder transforms the input waveform into a high-dimensional feature representation that is optimized to separate different sources. Then, the separation step computes a multiplicative mask for each of the target sources. Finally, the decoder reconstructs the original waveforms using the masked features.

Universal sound separation

The second type of method we tested are the universal sound separation methods TDCN++ and NMF, which can separate arbitrary sounds, not just musical instruments. This method more closely matches our problem of separating non-musical target sounds.

Non-negative matrix factorization (NMF) [START_REF] Cichocki | Fast local algorithms for large scale nonnegative matrix and tensor factorizations[END_REF][START_REF] Févotte | Algorithms for nonnegative matrix factorization with the betadivergence[END_REF]) is an unsupervised learning method, which takes an input non-negative matrix, in this case a spectrogram S, in R N ×M and decomposes it into a basis matrix A in R N ×J and a component matrix X in R J×M such that S ≈ A × X. We use the coordinate descent algorithm from the scikit-learn Python package based on [START_REF] Cichocki | Fast local algorithms for large scale nonnegative matrix and tensor factorizations[END_REF] with J set to 4, which selects for 4 features in the dataset. We then output 4 mono waveforms containing each isolated feature F i , where

F i = A i × X i for i = 1, 2, 3, or 4.
TDCN++ is a model proposed by [START_REF] Kavalerov | Universal sound separation[END_REF], based on Conv-Tasnet. It was modified in order to be used for general sound separation, and not only for speech. The main architecture remains the same as in [START_REF] Luo | Tasnet: Surpassing ideal time-frequency masking for speech separation[END_REF], and the changes only affect the masking network. More specifically, feature-wise normalization replaces global normalization, residual connection ranges are extended to improve the flow of information, and learnable scale parameters are introduced to weight the outputs of each residual layer. Those modifications help the model to better extract features without drastically increasing its size. We use a trained version of this model, which was trained on a custom dataset made of samples from movie production recordings from the Pro Sound Effects Library database (Pro Sound Effects, 2021). We use the model given as a baseline for the FUSS dataset [START_REF] Wisdom | What's all the fuss about free universal sound separation data?[END_REF], which is able to separate mixtures with a variable number of sources, ignoring the outputs that have a power below a given threshold.

Testing

The procedure we use for testing is as follows:

1. The target is created as a combination of four randomly chosen sources, which are offset to begin playing at different times.

2. The target, without any separation performed, is orchestrated using the entire orchestra.

This creates what we call the full target orchestration.

3. The full orchestra is randomly split into four equal sized subsets, each containing 7 to 8 instruments.

4. For each separation method:

1. The separation method is applied to the target, splitting the target into four sub-targets.

2. The separation is evaluated using the cosine distance.

3. The four sub-targets are separately orchestrated, each using a different subset of the orchestra. 4. The four orchestrations are combined to play simultaneously, creating the separated orchestration. 5. The distance between the target and the separated orchestration is calculated.

The orchestrations performed are static orchestrations, meaning the orchestration does not change over time; a single onset of notes is created for each sub-target, no matter if the target itself has multiple onsets. The OrchideaSOL database of orchestral samples is used with Orchidea to create the orchestrations [START_REF] Cella | OrchideaSOL: a dataset of extended instrumental techniques for computer-aided orchestration[END_REF]. This database contains recordings of extended playing techniques, which better fits the often noisy or inharmonic nature of our targets. The full orchestra, of which non-overlapping subsets were selected to orchestrate the sub-targets, contains the following instruments: 8 violins, 4 violas, 3 cellos, 1 bass, 2 oboes, 2 flutes, 2 clarinets in B♭, 2 bassoons, 2 trumpets, 2 trombones, and 2 French horns.

Evaluation

We compare the effectiveness of the separation methods in two ways. First, we evaluate the separation accuracy by calculating the cosine distance between the target and the output of a separation method. Second, we compare how each separation method affects the quality of the orchestration by calculating the spectral distance between the target and the separated orchestration. We use these distance metrics as they are specific to our problem of computer-assisted orchestration [START_REF] Abreu | Computer-aided musical orchestration using an artificial immune system[END_REF].

Separation evaluation

The source separation step of our experiments was evaluated independently from the orchestration process. To evaluate the source separation quality of an estimate-reference pair, we first frame the waveforms into non-overlapping segments of approximately 90ms. Then, for each frame, we compute the normalized magnitude spectrum and compare them using the cosine distance. Therefore, the distance between an estimated sub-target x and a reference x is given by Eq. 1.

D(x, x) = 1 n n i (1 - ⟨x i , xi ⟩ ∥x i ∥ • ∥ xi ∥ ) (1)
where x i is the normalized magnitude spectrum of the i-th frame. To align the estimates from the source separation models with their corresponding reference sub-targets, we consider all referencesestimations permutations, and we keep the best result.

Orchestration evaluation

Our goal is to evaluate the effect of source separation methods on orchestration. We determine the quality of an orchestration by calculating the spectral distance between the target and orchestration.

The spectral distance metric is proposed in Cella (2020) as part of the cost function used in Orchidea during the optimization, and used in Cella et al. (2020a) to compare accuracies of orchestrations. The equation takes in the normalized magnitude spectrum of the target x and of the solution x. Then for each bin k of the spectrum, it calculates the absolute difference between the amplitudes. The differing values of λ 1 and λ 2 allow the metric to penalize the solution in different ways. For our purposes, we used λ 1 = 0.5 and λ 2 = 10, which penalizes a solution that overshoots the harmonic energy of the target.

d(x, x) = λ 1 k δ k1 (x k -xk ) + λ 2 k δ k2 |x k -xk | (2)
where δ k1 = 1 if x k ≥ xk , 0 otherwise; and δ k2 = 1 if x k < xk , 0 otherwise.

We found that sometimes one or more of the outputs of a separation method would be silence. This happens when the method cannot distinguish between two or more sources and incorrectly identifies them as the same source. Therefore, some sub-targets can have multiple sources and some sub-targets can be silence. However, Orchidea still attempts to orchestrate silence, so when a sub-target is identified as silence, the orchestration of that sub-target is not included in the separated orchestration.

We identified silent sub-targets by calculating their root mean square (RMS) value. If the RMS was below a threshold of 0.05 then the sub-target orchestration was not included in the separated orchestration.

RM S(x) = 1 n n i x 2 i (3)
To determine whether source separation improves orchestration, we compare the distance between the target and the full target orchestration to the distance between the target and the separated orchestration. If the separated orchestration distance is lower, then performing separation improved the orchestration.

Results and discussion

The separation evaluation and orchestration evaluation results, averaged across 300 targets, are displayed in Tab. 1. The results show that performing source separation followed by a static orchestration of each individual sub-target improves the orchestration for all methods tested. Furthermore, when comparing the results between the separation evaluation and orchestration evaluation, we find that as the quality of the separation increases, so does the quality of the orchestration. In fact, the relative rankings between each method are exactly the same when evaluating separation and orchestration. NMF shows the best results, it has the lowest cosine distance for the separation evaluation and shows a 39% decrease in spectral distance compared to the full target orchestration.

A surprising result is the difference in accuracy for the two universal source separation methods: NMF and TDCN++. For both the separation evaluation and the orchestration evaluation, NMF performed the best out of all methods and TDCN++ performed the worst. We expect universal source separation to fit our problem better than music source separation, since our targets contain many non-musical sounds. However, our data suggests that the aspect of the separation method that affects orchestration is not whether it is universal or music source separation, but whether it is a supervised or unsupervised method.

The one unsupervised method we tested, NMF, outperformed all the supervised methods. This suggests that either unsupervised methods are better for our task, or that the data the supervised methods were trained on does not fit our problem well. It is important to note that for the supervised methods, we used pre-trained versions of the neural models. All of these models, except for TDCN++, were trained on musical data. However, we tested them on a wide range of targets, many of which were not strictly musical. This could explain why Demucs, Open-Unmix, and Conv-Tasnet performed worse than NMF. Similarly, one reason that TDCN++ may have performed worse than the other supervised methods is because it was trained on a different database than the other three supervised methods.

Conclusions and future work

In this paper, we addressed the potential of using source separation as a pre-processing step for computer-assisted orchestration. Based on the concept of dovetailing, our hypothesis is that separating a multi-layered sound and orchestrating the individual layers should result in an orchestration that better matches the timbre of the target. Our data confirms our hypothesis: all of the separation methods tested improve the resulting orchestration. Furthermore, there is a correlation between the effectiveness of the separation and the quality of the orchestration: as separation improves, so does the orchestration. This shows that our approach and metrics are consistent.

We also performed a qualitative evaluation of the orchestrations. An interesting result that we discovered from acoustic inspection is that the orchestrations generated after source separation tend to favor dovetailing. Since only a portion of the orchestra is assigned to each source and since the sources are inherently asynchronous with each other, the final orchestration is more fluid and appears more interesting from a musical standpoint.

More work is needed to identify why TDCN++, a universal source separation method, performed worse than the musical source separation models. Its performance is likely a result of the data it was trained on, which was not necessarily optimized for the sounds common to computer-assisted orchestration. A next step could be training this architecture ourselves on data that fits our problem better.

Another aspect that could be improved is the assignment of the orchestras. Currently, the subset of the orchestra that is assigned to each sub-target is randomly selected. However, it is possible that this random assignment could reduce the quality of the final orchestration. If, for example, one of the sub-targets contains mostly low pitch content, but the orchestra assigned to it contains mostly high-pitched instruments, the resulting orchestration will suffer. To solve this problem, a step could be added that jointly-optimizes the orchestras assigned to each sub-target, ultimately choosing an assignment that leads to the best orchestration of each sub-target.

As we stated in Section 1, we performed only static orchestration in our experiments in order to better observe the effects of source separation. However, dynamic orchestration is an important feature of computer-assisted orchestration and could greatly benefit from source separation. Future work could include applying source separation to dynamic targets and performing dynamic orchestration.

Finally, one or multiple of the separation methods we tested should be implemented in Orchidea in order to improve the orchestrations Orchidea is capable of creating.

Figure 1 :

 1 Figure 1: Diagram comparing the process of orchestration with and without separation. At the top, the full target orchestration is created simply by orchestrating the target. At the bottom, the separation process is performed, and the resulting sub-targets are individually orchestrated with subsets of the orchestra before being recombined to obtain the separated orchestration.

Figure 2 :

 2 Figure 2: Diagram showing how the full target orchestration is a single orchestration with the full orchestra, and the separated orchestration is four orchestrations with four subsets of the orchestra.

Table 1 :

 1 Results for separation and orchestration evaluation averaged across 300 targets. Separation evaluations uses the cosine distance defined in Eqn. 1. Orchestration evaluation uses the spectral distance defined in Eqn. 2. "Full target" corresponds to the orchestration of the target with no separation performed (seeSec. 3.3). Lower values are better for both metrics.

	Separation method Separation evaluation Orchestration evaluation
	TDCN++	0.924	8.59
	Conv-Tasnet	0.816	8.13
	Open-Unmix	0.668	7.21
	Demucs	0.655	6.90
	NMF	0.633	5.61
	Full target	-	9.24