Daniel Li

Hervé Queffélec

Luis Rodríguez-Piazza

Estimates for approximation numbers of some classes of composition operators on the Hardy space

Keywords: Mathematics Subject Classification 2010. Primary: 47B06 -Secondary: 30J10, 47B33 approximation numbers, Blaschke product, composition operator, cusp map, Hardy space, modulus of continuity, Schatten classes

We give estimates for the approximation numbers of composition operators on H 2 , in terms of some modulus of continuity. For symbols whose image is contained in a polygon, we get that these approximation numbers are dominated by e -c √ n . When the symbol is continuous on the closed unit disk and has a domain touching the boundary non-tangentially at a finite number of points, with a good behavior at the boundary around those points, we can improve this upper estimate. A lower estimate is given when this symbol has a good radial behavior at some point. As an application we get that, for the cusp map, the approximation numbers are equivalent, up to constants, to e -c n/ log n , very near to the minimal value e -c n . We also see the limitations of our methods. To finish, we improve a result of O. El-Fallah, K. Kellay, M. Shabankhah and H. Youssfi, in showing that for every compact set K of the unit circle T with Lebesgue measure 0, there exists a compact composition operator C ϕ : H 2 → H 2 , which is in all Schatten classes, and such that ϕ = 1 on K and |ϕ| < 1 outside K.

Introduction and notation

If the approximation numbers of some classes of operators on Hilbert spaces are well understood (for example, those of Hankel operators: see [START_REF] Megretskii | The inverse spectral problem for self-adjoint Hankel operators[END_REF]), it is not the case of those of composition operators. Though their behavior remains mysterious, some recent results are obtained in [START_REF] Li | On approximation numbers of composition operators[END_REF] and [START_REF] Lefèvre | Some new properties of composition operators associated with lens maps[END_REF] for approximation numbers of composition operators on the Hardy space H 2 . In [START_REF] Li | On approximation numbers of composition operators[END_REF], it is proved that one always has a n (C ϕ ) e -c n for some c > 0 ( [START_REF] Li | On approximation numbers of composition operators[END_REF], Theorem 3.1) and that this speed of decay can only be got when the symbol ϕ maps the unit disk D into a disk centered at 0 of radius strictly less than 1, i.e. ϕ ∞ < 1 ([14], Theorem 3.4).

In this paper, we give estimates which are somewhat general, in terms of some modulus of continuity. In Section 2, we obtain an upper estimate when the symbol ϕ is continuous on the closed unit disk and has an image touching non-tangentially the unit circle at a finite number of points, with a good behavior on the boundary around this point. As an application, we show that for symbols ϕ whose image is contained in a polygon a n (C ϕ ) ≤ a e -b √ n , for some constants a, b > 0; this has to be compared with [START_REF] Lefèvre | Some new properties of composition operators associated with lens maps[END_REF], Proposition 2.7, where it is shown that if ϕ is a univalent symbol such that ϕ(D) contains an angular sector centered on the unit circle and with opening θπ, 0 < θ < 1, then a n (C ϕ ) ≥ a e -b √ n , for some (other) positive constants a and b, depending only on θ. In Section 3, we obtain a lower bound when ϕ has a good radial behavior at the contact point. Both proofs use Blaschke products. This allows to recover the estimation a n (C λ θ ) ≈ e -c √ n obtained in [START_REF] Li | On approximation numbers of composition operators[END_REF], Proposition 6.3, and [START_REF] Lefèvre | Some new properties of composition operators associated with lens maps[END_REF], Theorem 2.1 for the lens map λ θ . In Section 4.1, we give another example, the cusp map, for which a n (C ϕ ) ≈ e -c n/ log n , very near the minimum value e -c n . We end that section by considering a one-parameter class of symbols, first studied by J. Shapiro and P. D. Taylor [START_REF] Shapiro | Compact, nuclear, and Hilbert-Schmidt composition operators on H 2[END_REF] and seeing the limitations of our methods. In Section 5, we improve a result of E.A. Gallardo-Gutiérrez and M.J. González (previously generalized by O. El-Fallah, K. Kellay, M. Shabankhah and H. Youssfi [START_REF] El-Fallah | Level sets and composition operators on the Dirichlet space[END_REF], Theorem 3.1). It is known that for every compact composition operator C ϕ : H 2 → H 2 , the set E ϕ = {e iθ ; |ϕ * (e iθ )| = 1} has Lebesgue measure 0. These authors showed ( [START_REF] Gallardo-Gutiérrez | Hausdorff measures, capacities and compact composition operators[END_REF]), with a rather difficult construction, that there exists a compact composition operator C ϕ : H 2 → H 2 such that the Hausdorff dimension of E ϕ is equal to 1 (and in [START_REF] El-Fallah | Level sets and composition operators on the Dirichlet space[END_REF], it is shown that for any negligible compact set K, there is a Hilbert-Schmidt operator C ϕ such that E ϕ = K). We improve this result in showing that for every compact set K of the unit circle T with Lebesgue measure 0, there exists a compact composition operator C ϕ : H 2 → H 2 , which is even in all Schatten classes, and such that E ϕ = K.

C ϕ : H 2 → H 2 by C ϕ (f ) = f • ϕ, called the composition operator of symbol ϕ.
Recall that if T : E → F is a bounded operator between two Banach spaces, the approximation numbers a n (T ) of T are defined by:

a n (T ) = inf{ T -R ; rank (R) < n} , n = 1, 2, . . .
The sequence a n (T ) n is non-increasing and, when F has the Approximation Property, T is compact if and only if a n (T ) tends to 0.

Definition 1.1 A modulus of continuity ω is a continuous function ω : [0, A] → R + ,
which is increasing, sub-additive, and vanishes at zero. Some examples are:

ω(h) = h α , 0 < α ≤ 1 ; ω(h) = h log 1 h ; ω(h) = 1 log 1 h •
For any modulus of continuity ω, there is a concave modulus of continuity ω ′ such that ω ≤ ω ′ ≤ 2ω (see [START_REF] Medvedev | On a concave differentiable majorant of a modulus of continuity[END_REF] for example); therefore we may and shall assume that ω is concave on [0, A]. In that case, ω -1 is convex, and

(1.1) r ω (x) := ω -1 (x) x is non-decreasing.
The notation u(t) v(t) means that u(t) ≤ Av(t) for some constant A > 0 and u(t) ≈ v(t) means that both u(t) v(t) and v(t) u(t).

Upper bound and boundary behavior

Definition 2.1 Let ω be a modulus of continuity and ϕ a symbol in the disk algebra A(D). Let ξ 0 ∈ ∂D ∩ ϕ(D). We say that the symbol ϕ has an ω-regular behavior at ξ 0 if, setting:

(2.1) γ(t) = ϕ(e it ) ,
and E ξ0 = {t ; γ(t) = ξ 0 }, there exists r 0 > 0 such that:

1) for some positive constant C > 0, one has, for every t 0 ∈ E ξ0 and |t-t 0 | ≤ r 0 :

(2.2) |γ(t) -γ(t 0 )| ≤ C 1 -|γ(t)| .
2) for some positive constant c > 0, one has, for for every t 0 ∈ E ξ0 and |t -t 0 | ≤ r 0 :

(2.3) c ω(|t -t 0 |) ≤ |γ(t) -γ(t 0 )| .
The first condition implies that the image of ϕ touches ∂D at the point ξ 0 , and non-tangentially. The second one implies that ϕ does not stay long near ξ 0 = γ(t 0 ).

Note that, due to (2.3), the intervals [t -r 0 /2, t + r 0 /2], for t ∈ E ξ0 are pairwise disjoint and therefore the set E ξ0 must be finite.

We shall make the following assumption (to avoid the Lipschitz class):

(2.4)

lim h→0 + ω(h) h = ∞ ; equivalently lim h→0 + ω -1 (h) h = 0.
Indeed, assume that γ is K-Lipschitz at some point

t 0 ∈ [0, 2π], namely |ϕ(e it ) -ϕ(e it0 )| ≤ K |t -t 0 |, with |ϕ(e it0 )| = 1; then m({t ∈ [0, 2π] ; |ϕ(e it ) -ϕ(e it0 )| ≤ h}) ≥ m({t ∈ [0, 2π] ; |t -t 0 | ≤ h/K}) = h/2πK ;
hence this measure in not o (h) and the composition operator C ϕ is not compact ( [START_REF] Maccluer | Compact composition operators on H p (B N )[END_REF], or [START_REF] Cowen | Composition Operators on Spaces of Analytic Functions[END_REF], Theorem 3.12).

In order to treat the case where the image of ϕ is a polygon, we need to generalize the above definition. We ask not only that ϕ is ω-regular at the points ξ 1 , . . . , ξ p of contact of ϕ(D) with ∂D, but a little bit more. Definition 2.2 Assume that ϕ(D) ∩ ∂D = {ξ 1 , . . . , ξ p }. We say that ϕ is globally-regular if there exists a modulus of continuity ω such that, writing E ξj = {t ; γ(t) = ξ j }, one has, for some r 1 , . . . , r p > 0

T = p j=1 E ξj + [-r j , r j ]
and for some positive constants C, c > 0, 1') one has, for j = 1, . . . , p, every t j ∈ E ξj and |t -t j | ≤ r j :

(2.5)

|γ(t) -γ(t j )| ≤ C 1 -|γ(t)| .
2') one has, for j = 1, . . . , p, every t j ∈ E ξj and |t -t j | ≤ r j :

(2.6) c ω(|t -t j |) ≤ |γ(t) -γ(t j )| .
Let us note that condition 1') is equivalent to say that ϕ(D) is contained in a polygon inside D whose vertices contain ξ 1 , . . . , ξ p , and these are the only vertices in the boundary ∂D. Of course, we may assume that (2.5) and (2.6) hold only when t is in a neighborhood of t j , since they will then hold for |t -t j | ≤ r j , provided we change the constants C, c. Before stating our theorem, let us introduce a notation. If ϕ is as in Definition 2.2 and σ, κ > 0 are some constants, we set:

(2.7) d N = σ log κ 2 -N ω -1 (κ 2 -N ) + 1,
where [ ] stands for the integer part. For every integer q ≥ 1, we denote by (2.8) N = N q the largest integer such that p N d N < q

(N q = 1 if no such N exists).
We then have the following result.

Theorem 2.3 Let ϕ be a symbol in A(D) whose image touches ∂D at the points ξ 1 , . . . , ξ p , and nowhere else. Assume that ϕ is globally-regular. Then, there are constants κ, K, L > 0, depending only on ϕ, such that, using the notation (2.7) and (2.8), one has, for every q ≥ 1:

(2.9)

a q (C ϕ ) ≤ K ω -1 (κ 2 -Nq ) κ 2 -Nq .
Before proving this theorem, let us indicate two applications. In these examples, we can give an upper estimate for all approximation numbers a n (C ϕ ), n ≥ 1 because we can interpolate between the integers N d N and (N + 1) d N +1 , which is not the case in general.

1) ω(h) = h θ , 0 < θ < 1, as this is the case for inscribed polygons (see the proof of the foregoing Theorem 2.4; here θ = max{θ 1 , . . . , θ p }, where θ 1 π, . . . , θ p π are the values of the angles of the polygon), as well as, with p = 2, for lens maps λ θ (see [START_REF] Shapiro | Composition Operators and Classical Function Theory, Universitext[END_REF], page 27, for the definition; see also [START_REF] Lefèvre | Some new properties of composition operators associated with lens maps[END_REF]). We have here

ω -1 (h) = h 1/θ . Hence d N ≈ N , N q ≈
√ q, and we then get from (2.9) that a q (C ϕ ) ≤ α 2 -δN for q N 2 , with δ > 0. Equivalently, for suitable constants α, β > 0,

(2.10) a n (C ϕ ) ≤ α e -β √ n , which is the result obtained in [12], Theorem 2.1. 2) ω(h) = 1 (log 1/h) α , 0 < α ≤ 1,
as this is the case, when α = 1, for the cusp map, defined below in Section 4.1 (with p = 1). Then, we have ω -1 (h) = e -h -1/α and d N ≈ 2 N/α , so that N q ≈ log q and 2 Nq/α ≈ q/ log q. Now, a simple computation gives:

(2.11) a n (C ϕ ) ≤ α e -βn/ log n .
Without assuming some regularity, one has the following general upper estimate.

Theorem 2.4 Let ϕ : D → D be an analytic self-map whose image is contained in a polygon P with vertices on the unit circle. Then, there exist constants α, β > 0, β depending only on P, such that:

(2.12) a n (C ϕ ) ≤ α e -β √ n .
In [START_REF] Lefèvre | Some new properties of composition operators associated with lens maps[END_REF], Proposition 2.7, it is shown that if ϕ is a univalent symbol such that ϕ(D) contains an angular sector centered on the unit circle and with opening θπ,

0 < θ < 1, then a n (C ϕ ) ≥ α e -β √
n , for some (other) positive constants α and β, depending only on θ. Note that the injectivity of the symbol is there necessary, since there exists (see the proof of Corollary 5.4 in [START_REF] Li | On approximation numbers of composition operators[END_REF]), for every sequence (ε n ) of positive numbers tending to 0, a symbol ϕ whose image is D \ {0}, and hence contains polygons), which is 2-valent, and for which a n (C ϕ ) e -εnn . This bound may be much smaller than e -β √ n .

Proof of Theorem 2.3. It follows the lines of that of [START_REF] Lefèvre | Some new properties of composition operators associated with lens maps[END_REF], Theorem 2.1.

Recall ([12], Lemma 2.4) that for every Blaschke product B with less than N zeros (each of them being counted with its multiplicity), one has:

(2.13) a N (C ϕ ) 2 sup 0<h<1,|ξ|=1 1 h S(ξ,h) |B(z)| 2 dm ϕ (z) ,
where S(ξ, h) = {z ∈ D ; |z -ξ| ≤ h} and m ϕ is the pull-back measure by ϕ of the normalized Lebesgue measure m on T. The proof will come from an adequate choice of a Blaschke product. Fix a positive integer N . Set, for j = 1, . . . , p and k = 1, 2, . . .:

(2.14)

p j,k = (1 -2 -k )ξ j
and consider the Blaschke product of length pN d (d being a positive integer, to be specified later) given by:

(2.15)

B(z) = p j=1 N k=1 z -p j,k 1 -p j,k z d .
Recall that we have set

(2.16) γ(t) = ϕ(e it ).
To use (2.13), note that if |γ(t) -ξ| ≤ h, then, for some j = 1, . . . , p and some t j ∈ E ξj , one has |t -t j | ≤ r j and, by (2.5), |γ(t) -

ξ j | ≤ C(1 -|γ(t)|) ≤ C |γ(t) -ξ| ≤ Ch.
Therefore, denoting by L j the number of elements of E ξj (which is finite by the remark following Definition 2.1):

a N (C ϕ ) 2 sup 0<h<1 1 h p j=1 L j {|γ(t)-ξj|≤Ch}∩{|t-tj|≤rj } |B[γ(t)]| 2 dt 2π ,
and we only need to majorize the integrals:

I j (h) = {|γ(t)-ξj|≤Ch}∩{|t-tj|≤rj} |B γ(t) | 2 dt 2π •
Moreover, it suffices, by interpolation, to do that with h = h n , where

h n = 2 -n . By (2.6), for |t -t j | ≤ r j and |γ(t) -ξ j | ≤ Ch n , one has c ω(|t -t j |) ≤ |γ(t) -ξ j | ≤ Ch n = C 2 -n , which implies that (2.17) |t -t j | ≤ ω -1 (c -1 C 2 -n ).

Let

(2.18)

s n = ω -1 (c -1 C 2 -n ).
One has:

I j (h n ) ≤ {|t-tj |≤sn}∩{|t-tj|≤rj } |B γ(t) | 2 dt 2π •
For n ≥ N , we simply majorize |B γ(t) | by 1 and we get:

1 h n I j (h n ) ≤ 1 h n 2s n 2π = c -1 C π 1 c -1 C 2 -n ω -1 (c -1 C 2 -n ) ≤ c -1 C π ω -1 (c -1 C 2 -N ) c -1 C 2 -N , since the function ω -1 (x)/x is non-decreasing.
When n ≤ N -1, we write:

I j (h n ) ≤ {|t-tj |≤sN }∩{|t-tj |≤rj} |B γ(t) | 2 dt 2π + {sN <|t-tj|≤sn}∩{|t-tj |≤rj} |B γ(t) | 2 dt 2π •
The first integral is estimated as above. For the second one, we claim that:

Claim 2.5 For some constant χ < 1, one has, for j = 1, . . . , p and every t j ∈ E ξj :

(2.19) |B γ(t) | ≤ χ d when |t -t j | > s N and |t -t j | ≤ r j .
To see that, we shall use [START_REF] Lefèvre | Some new properties of composition operators associated with lens maps[END_REF], Lemma 2.3. Let us recall that this lemma asserts that for w, w 0 ∈ D satisfying |w -w 0 | ≤ M min(1 -|w|, 1 -|w 0 |) for some positive constant M , one has:

(2.20) w -w 0 1 -w 0 w ≤ M √ M 2 + 1 • Let t such that |t -t j | ≤ r j and |t -t j | > s N . We have, on the one hand, ω(|t -t j |) ≥ ω(s N ) = c -1 C 2 -N ,
and, on the other hand, since |γ(t

j )| = |ξ j | = 1 c ω(|t -t j |) ≤ |γ(t) -γ(t j )| ≤ C(1 -|γ(t)|) ; hence 1 -|γ(t)| ≥ 2 -N . Let 1 ≤ k ≤ N such that 2 -k ≤ 1 -|γ(t)| < 2 -k+1 . Since |p j,k | = 1 -2 -k , we have: |γ(t) -p j,k | ≤ |γ(t) -ξ j | + |ξ j -p j,k | ≤ C(1 -|γ(t)|) + 2 -k ≤ (2C + 1)2 -k . Hence |γ(t) -p j,k | ≤ M min 1 -|γ(t)|, 1 -|p j,k | , with M = 2C + 1. By (2.20), we get γ(t)-p jk 1-p j,k γ(t) ≤ χ, where χ = M/ √ M 2 + 1 is < 1, and therefore |B[γ(t)]| ≤ χ d .
We can now end the proof of Theorem 2.3. We get:

1 h n {sN <|t-tj|≤sn}∩{|t-tj |≤rj} |B γ(t) | 2 dt 2π ≤ 1 h n 2s n 2π χ 2d = 1 h n ω -1 (c -1 C 2 -n ) π χ 2d = c -1 C π ω -1 (c -1 C 2 -n ) c -1 C 2 -n χ 2d ≤ 1 π ω -1 (c -1 C) χ 2d , since ω -1 (x)/x is non-decreasing. We therefore get, setting κ = c -1 C and L = L 1 + • • • + L p : 1 h n p j=1 L j {|γ(t)-ξj|≤Chn}∩{|t-tj |≤rj} |B[γ(t)]| 2 dt 2π ≤ κL π ω -1 (κ 2 -N ) κ 2 -N + L ω -1 (κ) π χ 2d .
Choose now d = d N , where d N is defined by (2.7), with σ = 1/ log(χ -2 ). Then χ 2d ≤ ω -1 (κ 2 -N )/(κ 2 -N ), and, since the Blaschke product B has now p N d N zeroes, we get, for some positive constant K:

a pN dN +1 (C ϕ ) ≤ K ω -1 (κ 2 -N ) κ 2 -N ,
and that ends the proof of Theorem 2.3.

Proof of Theorem 2.4. It suffices to consider the case when ϕ is a conformal map from D onto P. Indeed, let ψ be such a conformal map. In the general case, our assumption allows to write ϕ = ψ • u, where

u = ψ -1 • ϕ : D → D is analytic. It follows that C ϕ = C u • C ψ and that a n (C ϕ ) ≤ C u a n (C ψ ).
Therefore, we may and shall assume that ϕ itself is this conformal map. Let us denote by ξ 1 , . . . , ξ p the vertices of P. Let 0 < πµ j < π be the exterior angle of P at ξ j , namely the complement to π of the interior angle; so that: p j=1 µ j = 2 , and 0 < µ j < 1.

If one sets θ j = 1 -µ j , one has 0 < θ j < 1.

We then use the explicit form of ϕ given by the Schwarz-Christoffel formula ( [START_REF] Nehari | Conformal mapping[END_REF], page 193):

(2.21) ϕ(z) = A z 0 dw (a 1 -w) µ1 • • • (a p -w) µp + B ,
for some constants A = 0 and B ∈ C and where a 1 , . . . , a p ∈ ∂D are such that ξ j = ϕ(a j ), j = 1, . . . , p. If, as before, we write γ(t) = ϕ(e it ), we have ξ j = γ(t j ), with a j = e itj (note that here E ξj = {t j }).

As we already said, condition (2.5) is trivially satisfied for a polygon.

To end the proof, we use Theorem 2.3 and its Example 1. For that it suffices to show that, for |t -t j | small enough, we have:

(2.22) |γ(t) -ξ j | ≈ |t -t j | θj .
If z ∈ D is close to a j , it follows from (2.21) that we can write

ϕ(z) = A z 0 f j (w) dw (a j -w) µj + B,
where f j is holomorphic near a j and f j (a j ) = 0 since

|f j (a j )| = k =j,1≤k≤p |a j -a k | -µ k .
Write f j (w) = f j (a j ) + (a j -w)g j (w) where g j is holomorphic near a j . We get:

ϕ(z) = Af j (a j ) z 0 dw (a j -w) µj + B + z 0 g j (w)(a j -w) θj dw := Af j (a j ) z 0 dw (a j -w) µj + B + ψ j (z),
which can still be written (since θ j > 0):

(2.23) ϕ(z) = λ j (a j -z) θj + c j + ψ j (z),
where λ j = 0, c j ∈ C, ψ j is Lipschitz near a j and ξ j = ϕ(a j ) = c j + ψ j (a j ). Now, we easily get (2.22). Indeed, for t near t j , it follows from (2.23) that (recall that γ(t) = ϕ(e it ) and γ(t j ) = ξ j ):

|γ(t) -γ(t j )| = |λ j | |e it -e itj | θj + O (|t -t j |),
which the claimed estimate (2.22) since λ j = 0 and |t-t j | is negligible compared to |t -t j | θj ≈ |e it -e itj | θj .

Lower bound and radial behavior

We shall consider symbols ϕ taking real values in the real axis (i.e. its Taylor series has real coefficients) and such that lim r→1 -ϕ(r) = 1, with a given speed. Definition 3.1 We say that the analytic map ϕ : D → D is real if it takes real values on ] -1, 1[, and that ϕ is an ω-radial symbol if it is real and there is a modulus of continuity ω : [0, 1] → [0, 2] such that:

(3.1) 1 -ϕ(r) ≤ ω(1 -r) , 0 ≤ r < 1 .
With those definitions and notations, one has: Theorem 3.2 Let ϕ be a real and ω-radial symbol. Then, for the approximation numbers a n (C ϕ ) of the composition operator C ϕ of symbol ϕ, one has the following lower bound:

(3.2) a n (C ϕ ) ≥ c sup 0<σ<1 ω -1 (a σ n ) a σ n exp - 20 1 -σ ,
where a = 1 -ϕ(0) > 0 and c is another constant depending only on ϕ.

Observe that, for the lens map λ θ (see [START_REF] Lefèvre | Some new properties of composition operators associated with lens maps[END_REF], Lemma 2.5), we have ω -1 (h) ≈ h 1/θ , so that adjusting σ = 1 -1/ √ n, we get

(3.3) a n (C λ θ ) ≥ c exp -C √ n ,
which is the result of [START_REF] Li | On approximation numbers of composition operators[END_REF], Proposition 6.3.

For the cusp map ϕ (see Section 4.1), we have ω -1 (h) ≈ e -C ′ /h , so that taking σ = exp(-log n/2n), we get:

(3.4) a n (C ϕ ) ≥ c exp(-C n/ log n).
We shall use the same methods as for lens maps (see [START_REF] Li | On approximation numbers of composition operators[END_REF], Proposition 6.3). We need a lemma. Recall (see [START_REF] Hoffman | Banach Spaces of Analytic Functions[END_REF] pages 194-195, or [START_REF] Nikolski | Operators, Functions and Systems: An Easy Reading[END_REF] pages 302-303) that if (z j ) is a Blaschke sequence, its Carleson constant δ is defined as δ = inf j≥1 (1 -|z j | 2 ) |B ′ (z j )|, where B is the Blaschke product whose zeros are the z j 's. Now (see [START_REF] Garnett | Bounded Analytic Functions, Revised first version[END_REF], Chapter VII, Theorem 1.1), every H ∞ -interpolation sequence (z j ) is a Blaschke sequence and its Carleson constant δ is connected to its interpolation constant C by the inequalities

(3.5) 1/δ ≤ C ≤ κ/δ 2
where κ is an absolute constant (actually

C ≤ κ 1 (1/δ)(1+log 1/δ)). Now, if (z j )
is a H ∞ -interpolation sequence with constant C, the sequence of the normalized reproducing kernels

f j = K zj / K zj satisfies C -1 |λ j | 2 1/2 ≤ λ j f j H 2 ≤ C |λ j | 2 1/2
(see [START_REF] Li | On approximation numbers of composition operators[END_REF], Lemma 2.2).

Lemma 3.3 Let ϕ : D → D be an analytic self-map. Let u = (u 1 , . . . , u n ) be a finite sequence in D and set v j = ϕ(u j ), v = (v 1 , . . . , v n ). Denote by δ v the Carleson constant of the finite sequence v and set

µ 2 n = inf 1≤j≤n 1 -|u j | 2 1 -|ϕ(u j )| 2 •
Then, for some constant c ′ > 0, we have the lower bound:

(3.6) a n (C ϕ ) ≥ c ′ δ 4 v µ n .
Proof. Recall first that the Carleson constant δ of a Blaschke sequence (z j ) is also equal to:

δ = inf k≥1 j =k ρ(z k , z j ) ,
where ρ(z, ζ) = z-ζ 1-z ζ is the pseudo-hyperbolic distance between z and ζ. Now, the Schwarz-Pick Lemma (see [START_REF] Beardon | The hyperbolic metric and geometric function theory, Quasiconformal mappings and their applications[END_REF], Theorem 3.2) asserts that every analytic selfmap of D contracts the pseudo-hyperbolic distance. Hence ρ ϕ(u j ), ϕ(u k ) ≤ ρ(u j , u k ) and so, if δ u and δ v denote the Carleson constants of u and v:

δ u ≥ δ v .
Let now R be an operator of rank < n. There exists a function f = n j=1 λ j K uj ∈ H 2 ∩ ker R with f = 1. We thus have:

C * ϕ -R 2 ≥ C * ϕ (f ) -R(f ) 2 2 = C * ϕ (f ) 2 2 = n j=1 λ j K vj 2 2 ≥ C -2 v n j=1 |λ j | 2 K vj 2 2 = C -2 v n j=1 |λ j | 2 1 -|v j | 2 ≥ C -2 v µ 2 n n j=1 |λ j | 2 1 -|u j | 2 ≥ C -2 u C -2 v µ 2 n f 2 2 = C -2 u C -2 v µ 2 n ≥ κ -4 δ 4 u δ 4 v µ 2 n ≥ κ -4 δ 8 v µ 2 n ,
and hence

a n (C ϕ ) ≥ κ -2 δ 4 v µ n .
Remark. This lemma allows to give, in the Hardy case, a simpler proof of Theorem 4.1 in [START_REF] Li | On approximation numbers of composition operators[END_REF], avoiding the use of Lemma 2.3 and Lemma 2.4 (concerning the backward shift) in that paper. Recall that this theorem says that for every non-increasing sequence (ε n ) n≥1 of positive real numbers tending to 0, there exists a univalent symbol ϕ such that ϕ(0) = 0 and C ϕ :

H 2 → H 2 is compact, but a n (C ϕ ) ε n for every n ≥ 1.
Let us sketch briefly the argument. We use the notation of [START_REF] Li | On approximation numbers of composition operators[END_REF], Lemma 4.6. The symbol ϕ is defined as ϕ(z) = σ -1 (e -1 σ(z)), where σ is some conformal map σ : D → Ω. We set

A j = (1/C 0 ) log(1/ε j+1 ), r j = σ -1 (e j ).
Then ϕ(r j+1 ) = r j and (see [START_REF] Li | On approximation numbers of composition operators[END_REF], pages 444-446):

1 -r j+1 1 -r j ≥ exp(-2C 0 A j ).
We shall apply the above Lemma 3.3 with u j = r j . Then v j = ϕ(u j ) = r j-1 . Hence

1 -|u j | 2 1 -|v j | 2 ≥ 1 2 1 -|u j | 1 -|v j | = 1 2 1 -r j 1 -r j-1 ≥ 1 2 exp(-2C 0 A j-1 ) = 1 2 ε 2 j ≥ 1 2 ε 2 n .
It follows that µ n ≥ ε n / √ 2. On the other hand, (r j ) j≥1 is an interpolating sequence (see [START_REF] Li | On approximation numbers of composition operators[END_REF], Lemma 4.6); hence there is a constant δ > 0 (which does not depend on n ≥ 1) such that δ v ≥ δ. Therefore Lemma 3.3 gives

a n (C ϕ ) ≥ c δ 4 ε n ,
which gives Theorem 4.1 of [START_REF] Li | On approximation numbers of composition operators[END_REF].

Proof of Theorem 3.2. Fix 0 < σ < 1 and define inductively u j ∈ [0, 1) by u 0 = 0 and the relation

1 -ϕ(u j+1 ) = σ[1 -ϕ(u j )] with 1 > u j+1 > u j
(using the intermediate value theorem).

Setting v j = ϕ(u j ), we have -1 < v j < 1,

(3.7) 1 -v j+1 1 -v j = σ, and
(3.8) 1 -v n = a σ n , with a = 1 -ϕ(0).
Now observe that, for 1 ≤ j ≤ n, one has, due to the positivity of u j and v j , to (3.1), and the fact that r ω (x) = ω -1 (x)/x is increasing:

1 -|u j | 2 1 -|v j | 2 ≥ 1 -u j 2(1 -v j ) ≥ 1 2 ω -1 (1 -v j ) 1 -v j = 1 2 r ω (1 -v j ) ≥ 1 2 r ω (1 -v n ) = 1 2 r ω (a σ n ),
which proves that µ 2 n ≥ r ω (a σ n )/2. Furthermore, the sequence (v j ) satisfies, by (3.7), a condition very similar to Newman's condition with parameter σ. In fact, for k > j, we have

|v k -v j | |1 -v k v j | = (1 -v j ) -(1 -v k ) (1 -v j ) + v j (1 -v k ) ≥ (1 -v j ) -(1 -v k ) (1 -v j ) + (1 -v k ) = 1 -σ k-j 1 + σ k-j •
Analogously, for j > k, we have

|v k -vj | |1-v k vj | ≥ 1-σ j-k 1+σ j-k .
Thus, as in the proof of [START_REF] Duren | Theory of H p -spaces[END_REF], Theorem 9.2, we have, for every k, [START_REF] Li | On approximation numbers of composition operators[END_REF], Lemma 6.4. Finally, use (3.6) to get:

j =k ρ(v j , v k ) = j =k |v k -v j | |1 -v k v j | ≥ ∞ l=1 1 -σ l 1 + σ l 2 . Consequently, δ v ≥ ∞ l=1 1-σ l 1+σ l 2 ≥ exp -5 1-σ , by
a n (C ϕ ) ≥ c ′ δ 4 v µ n ≥ c exp - 20 1 -σ r ω (a σ n ).
Taking the supremum over σ, that ends the proof of Theorem 3.2.

Remark. The proof shows that (3.9) a n (C ϕ ) ≥ sup u1,...,un∈(0,1)

inf f ∈ Ku 1 ,...,Ku n f =1 C * ϕ f ,
where K u1 , . . . , K un is the linear space generated by n distinct reproducing kernels K u1 , . . . , K un . But if B is the Blaschke product with zeros u 1 , . . . , u n , then K u1 , . . . , K un = (BH 2 ) ⊥ , the model space associated to B. Hence

(3.10) a n (C ϕ ) ≥ sup B inf f ∈(BH 2 ) ⊥ f =1 C * ϕ f ,
where the supremum is taken over all Blaschke products with n zeros on the real axis (0, 1). This has to be compared with the upper bound (which gives (2.13), see [START_REF] Lefèvre | Some new properties of composition operators associated with lens maps[END_REF], proof of Lemma 2.4):

(3.11) a n (C ϕ ) ≤ inf B C ϕ |BH 2 = inf B sup f ∈BH 2 f =1 C ϕ f ,
where the infimum is over the Blaschke products with less than n zeros (in the Hilbert space H 2 , the approximation number a n (C ϕ ) is equal to the Gelfand number c n (C ϕ ), which is, by definition, less or equal to This map was first introduced in [START_REF] Lefèvre | Compact composition operators on Bergman-Orlicz spaces[END_REF] (see also [START_REF] Li | Compact composition operators on Hardy-Orlicz and Bergman-Orlicz spaces[END_REF]). Explicitly, ϕ is defined as follows.

C ϕ |BH 2 , since BH 2 is of codimension < n).

Examples

The cusp map

We first map D onto the half-disk D + = {z ∈ D ; Re z > 0}. To do that, map D onto itself by z → iz; then map D onto the upper half-plane H = {z ∈ C ; Im z > 0} by:

T (u) = i 1 + u 1 -u • Take the square root to map H in the first quadrant Q 1 = {z ∈ H ; Re z > 0},
and go back to the half-disk {z ∈ D ; Im z < 0} by T -1 : T -1 (s) = 1+is is-1 ; finally, make a rotation by i to go onto D + . We get:

(4.1) ϕ 0 (z) = z -i iz -1 1/2 -i -i z -i iz -1 1/2 + 1 • One has ϕ 0 (1) = 0, ϕ 0 (-1) = 1, ϕ 0 (i) = -i and ϕ 0 (-i) = i. The half-circle {z ∈ T ; Re z ≥ 0} is mapped onto the segment [-i, i] and the segment [-1, 1] onto the segment [0, 1]. Set now, successively, (4.2)
ϕ 1 (z) = log ϕ 0 (z), ϕ 2 (z) = - 2 π ϕ 1 (z) + 1, ϕ 3 (z) = 1 ϕ 2 (z) ,
and finally:

(4.3) ϕ(z) = 1 -ϕ 3 (z) .
Hence:

(4.4) 1 -ϕ(z) = 1 1 + 2 π log 1/|ϕ 0 (z)| -i 2 π arg ϕ 0 (z) • ϕ 2 maps D onto the semiband {z ∈ C ; Re z > 1 and |Im z| < 1}. One has ϕ(1) = 1, ϕ(-1) = 0, ϕ(i) = (1 + i)/2 and ϕ(-i) = (1 -i)/2.
The domain ϕ(D) is edged by three circular arcs of radii 1/2 and of respective centers 1/2, 1+i/2 and 1-i/2. The real interval ]-1, 1[ is mapped onto the real interval ] 0, 1[ and the half-circle {e iθ ; |θ| ≤ π/2} is sent onto the two circular arcs tangent at 1 to the real axis. Hence, when r tends to 1 -, one has:

(4.6) 1 -ϕ(r) ∼ π 2 1 log(1/γ) ∼ π 2 1 log(1/(1 -r)) •
2) For |θ| < π/2, one has:

(4.7) ϕ 0 (e iθ ) = -i tan(θ/2) 1 + 1 -tan 2 (θ/2) •
Hence, when θ tends to 0, one has:

(4.8) 1 -ϕ(e iθ ) ∼ π 2 1 log(1/|θ|) •
Proof. 1) One has:

T (ir) = r -i ir -1 = - 2r 1 + r 2 + i 1 -r 2 1 + r 2 = -sin α + i cos α ,
with r = tan(α/2); hence T (ir) = cos(α + π/2) + i sin(α + π/2) = e i(α+π/2) . Set β = α 2 + π 4 ; one gets:

ϕ 0 (r) = e iβ -i -ie iβ + 1 = cos β 1 + sin β = sin γ 1 + cos γ = tan(γ/2)
with γ = (π/2) -β = (π/4) -(α/2) = (π/4) -tan -1 r. Then (4.6) follows.

2) Let τ = π 2 -θ; one has:

T (ie iθ ) = e iθ -i ie iθ -1 = -cos θ 1 + sin θ = -sin τ 1 + cos τ = -tan(τ /2).
Note that 0 < τ /2 < π/2 since |θ| < π/2; hence tan(τ /2) > 0. Therefore:

ϕ 0 (e iθ ) = i tan(τ /2) -i -i.i tan(τ /2) + 1 = i tan(τ /2) -1 tan(τ /2) + 1 • But tan(τ /2) = tan π 4 - θ 2 = 1 -tan(θ/2) 1 + tan(θ/2) ;
it follows that:

ϕ 0 (e iθ ) = i 1 -tan(θ/2) -1 + tan(θ/2) 1 -tan(θ/2) + 1 + tan(θ/2) = i 1 -tan(θ/2) -1 + tan(θ/2) 1 -tan(θ/2) + 1 + tan(θ/2) 2 = -i tan(θ/2) 1 + 1 -tan 2 (θ/2) •
Now, since ϕ 0 (e iθ ) ∼ -iθ/4 as θ tends to 0, we get (4.8).

It follows from this lemma and from Theorem 2.3 and Theorem 3.2 that one has the following estimate.

Theorem 4.3 For the approximation numbers a n (C ϕ ) of the composition operator C ϕ : H 2 → H 2 of symbol the cusp map ϕ, we have: (4.9)

e -c1 n/ log n a n (C ϕ ) e -c2 n/ log n , n = 2, 3, . . . , for some constants c 1 > c 2 > 0.

Proof.

The Shapiro-Taylor map

This one-parameter map ς θ , θ > 0, was introduced by J. Shapiro and P. Taylor in 1973 ( [START_REF] Shapiro | Compact, nuclear, and Hilbert-Schmidt composition operators on H 2[END_REF]) and was further studied, with a slightly different definition, in [START_REF] Lefèvre | Some examples of compact composition operators on H 2[END_REF], Section 5. J. Shapiro and P. Taylor proved that C ς θ : H 2 → H 2 is always compact, but is Hilbert-Schmidt if and only if θ > 2. It is proved in [START_REF] Lefèvre | Some examples of compact composition operators on H 2[END_REF], Theorem 5.1, that C ς θ is in the Schatten class S p if and only if p > 4/θ.

Here, we shall use these maps ς θ to see the limitations of our previous methods.

We first recall their definition. For ε > 0, we set V ε = {z ∈ C ; Re z > 0 and |z| < ε}. For ε = ε θ > 0 small enough, one can define (4.10)

f θ (z) = z(-log z) θ ,
for z ∈ V ε , where log z will be the principal determination of the logarithm. Let now g θ be the conformal mapping from D onto V ε , which maps T = ∂D onto ∂V ε , defined by g θ (z) = ε ϕ 0 (z), where ϕ 0 is given in (4.1). Then, we define:

(4.11) ς θ = exp(-f θ • g θ ).
One has ς θ (1) = 1 and g θ (e it ) ∼ -it/4 as t tends to 0, by Lemma 4.2; hence, when t is near of 0:

|1 -ς θ (e it )| ≈ |f θ [g θ (e it )]| ≈ |t| [log(1/|t|)] θ .
If we were allowed to apply Theorem 2.3, we would get that a n (C ς θ ) 1/n θ/4 , which would be in accordance with the fact that C ς θ is in the Schatten class S p if and only if p > 4/θ. However, condition (2.2) is not satisfied: by [START_REF] Lefèvre | Some examples of compact composition operators on H 2[END_REF], equations (5.5) and (5.6), one has

1 -|ς θ (e it )| ≈ |t|(log 1/|t|) θ-1 , whereas |1 -ς θ (e it )| ≈ |t|(log 1/|t|) θ .
On the other hand, by the Lemma 4.2 again, g θ (r) ∼ ε(1 -r)/4 as r tends to 1; hence, when r is near to 1:

1 -ς θ (r) ≈ (1 -r) log 1/(1 -r) θ ,
so ς θ is a real ω-radial symbol with ω(t) = t(log 1/t) θ . Hence, we get from Theorem 3.2:

a n (C ς θ ) 1 n θ/2
, taking σ = 1/e in (3.2). However, this lower estimate is not the right one, since C ς θ is in S p if and only if p > 4/θ.

Contact points

It is well-known (and easy to prove) that for every compact composition operator C ϕ : H 2 → H 2 , the set of contact points

E ϕ = {e iθ ; |ϕ * (e iθ )| = 1}
has Lebesgue measure 0. A natural question is: to what extent is this negligible set arbitrary? The following partial answer was given by E.A. Gallardo-Gutiérrez and M.J. González in [START_REF] Gallardo-Gutiérrez | Hausdorff measures, capacities and compact composition operators[END_REF]. As an application of our previous results, we shall extend these results, with a very simple proof. Our composition operator will not even be compact, or Hilbert-Schmidt, but in all Schatten classes S p , and moreover its approximation numbers will be as small as possible.

Theorem 5.3 Let K be a Lebesgue-negligible compact set of the circle T. Then, there exists a Schur function ψ ∈ A(D), the disk algebra, such that E ψ = K, ψ(e iθ ) = 1 for all e iθ ∈ K, and: Actually, we can improve on the previous theorem by proving the following result. This result is optimal because if ψ ∞ = 1, we know (see [START_REF] Li | On approximation numbers of composition operators[END_REF], Theorem 3.4) that lim inf n→∞ [a n (C ψ )] 1/n = 1, so we cannot hope to get rid with the forthcoming vanishing sequence (ε n ) n .

Theorem 5.4 Let K be a Lebesgue-negligible compact set of the circle T and (ε n ) n a sequence of positive real numbers with limit zero. Then, there exists a Schur function ϕ ∈ A(D) such that E ϕ = K, ϕ(e iθ ) = 1 for all e iθ ∈ K, and Lemma 5.5 Let δ be a nondecreasing positive function on (0, 1] tending to 0 as h → 0. Then, there exists a Schur function ψ ∈ A(D) such that ψ(1) = 1, |ψ(ξ)| < 1 for ξ ∈ T \ {1}, and such that ρ ψ (h) ≤ δ(h), for h > 0 small enough.

Once we have the lemma, in view of the upper bound in [START_REF] Li | On approximation numbers of composition operators[END_REF], Theorem 5.1, for approximation numbers, we can adjust the function δ so as to have a n (C ψ ) ≤ Ke -nεn . Then, we compose ψ with a peaking function χ as in the previous section and the map ϕ = ψ • χ fulfills the requirements of Theorem 5.4, with C = K χ .

Notation.

  We denote by D the open unit disk and by T = ∂D the unit circle; m is the normalized Lebesgue measure on T: dm(t) = dt/2π. The disk algebra A(D) is the space of functions which are continuous on the closed unit disk D and analytic in the open unit disk. If H 2 is the usual Hardy space on D, every analytic self-map ϕ : D → D (also called Schur function) defines, by Littlewood's subordination principle, a bounded operator

Definition 4 . 1

 41 The cusp map is the conformal mapping ϕ sending the unit disk D onto the domain represented on Figure1.

Figure 1 :

 1 Figure 1: Cusp map domain

Lemma 4.2 1 )

 1 For 0 < r < 1, let γ = π 4 -arctan r = arctan[(1 -r)/(1 + r)]; then: (4.5) ϕ 0 (r) = tan(γ/2) .

Theorem 5 . 1 (

 51 E.A. Gallardo-Gutiérrez and M.J. González) There is a compact composition operator C ϕ on H 2 such that the Hausdorff dimension of E ϕ is one.This was generalized by O. El-Fallah, K. Kellay, M. Shabankhah, and H. Youssfi ([START_REF] El-Fallah | Level sets and composition operators on the Dirichlet space[END_REF], Theorem 3.1): Theorem 5.2 (O. El-Fallah, K. Kellay, M. Shabankhah, H. Youssfi) For every compact set K of measure 0 in T, there exists a Schur function ϕ ∈ A(D), the disk algebra, such that the associated composition operator C ϕ is Hilbert-Schmidt on H 2 and E ϕ = K.

( 5 . 1 )

 51 a n (C ψ ) ≤ a exp(-b n/ log n).In particular, C ψ ∈ p>0 S p .Proof. According to the Rudin-Carleson theorem ([START_REF] Bishop | A general Rudin-Carleson theorem[END_REF]), we can find χ ∈ A(D) such that χ = 1 on K and |χ| < 1 on D \ K.Consider now the cusp map ϕ, defined in Section 4.1. One has ϕ ∈ A(D), ϕ(1) = 1 anda n (C ϕ ) ≤ a ′ exp(-b n/ log n).We now spread the point 1 by composing with the function χ, which is equal to 1 on the whole of K. We check that the composed map ψ = ϕ • χ has the required properties.That ψ ∈ A(D) is clear. For z ∈ K, one has ψ(z) = ϕ(1) = 1, and for z ∈ D \ K, one has |χ(z)| < 1; hence |ψ(z)| < 1.To finish, sinceC ψ = C χ • C ϕ , have a n (C ψ ) ≤ C χ a n (C ϕ ) ≤ a ′ C χ exp(-b n/ log n) := σ n ,proving the result (with a = a ′ C χ ), since clearly ∞ n=1 σ p n < ∞ for each p > 0.

( 5 . 2 )

 52 a n (C ϕ ) ≤ C exp(-n ε n ) ,where C is a positive constant.This theorem is a straightforward consequence of the following lemma. Recall that the Carleson function of the Schur function ψ : D → D is defined by:ρ ψ (h) = sup|ξ|=1 m({t ∈ T ; |ψ(e it )| ≥ 1 -h and | arg(ψ(e it ) ξ)| ≤ πh}).

  1) Upper estimate. Note first that, since the domain ϕ(D) is contained in the right half-plane and in the symmetric angular sector of vertex 1 and opening π/2, there is a constant C > 0 such that |1 -γ(t)| ≤ C (1 -|γ(t)|) and we have (2.2). Then (4.8) in Lemma 4.2 gives (2.3). The upper estimate is hence given in Theorem 2.3 and (2.11). 2) Lower estimate. By Lemma 4.2, (4.6), one has (3.1). Since ϕ is a real symbol, the upper estimate follows from Theorem 3.2, and (3.4).

* Supported by a Spanish research project MTM 2009-08934.

Proof of Lemma 5.5. We use a slight modification of the map g constructed in [START_REF] Lefèvre | Some revisited results about composition operators on Hardy spaces[END_REF], pages 66-67. Instead of taking a conformal map from D to the domain used in [START_REF] Lefèvre | Some revisited results about composition operators on Hardy spaces[END_REF], we modify this domain by limiting it to the right-hand side (by, say, a semicircle), as on the Figure 2. Let Ω this domain. This domain is limited by the two hyperbolas y = 1/x and y = (1/x) + 4π. The limiting semicircle is chosen in order that Im w ≥ 1 for w ∈ Ω. The lower part of the "saw-teeth" have an imaginary part equal to 4πn. If a ∈ Ω is fixed and Ω n is the part of the domain Ω such that Im w < 4πn, the horizontal sizes of the "saw-teeth" are chosen in order that the harmonic measure ω Ω (a, ∂Ω ∂Ω n ) is ≤ δ n := δ(1/16π(n + 1)). Note that ∂Ω \ ∂Ω n ⊇ {w ∈ ∂Ω ; Im w > 4πn} (see [START_REF] Lefèvre | Some revisited results about composition operators on Hardy spaces[END_REF], Lemma 4.2).

Figure 2: Domain Ω By Carathéodory-Osgood's Theorem (see [START_REF] Palka | An Introduction to Complex Function Theory[END_REF], Theorem IX.4.9), there is a unique homeomorphism g from D onto Ω ∪ {∞} which maps conformally D onto Ω and such that g(0) = a and g(1) = ∞ (we may choose these two values because if h : D → Ω ∪ {∞} is such a map, and u is the automorphism of D such that u(0) = h -1 (a) and u(1) = h -1 (∞), then g = h • u suits -alternatively, having choosen h(0) = a, then, if h(e iθ0 ) = ∞, we take g(z) = h(e iθ0 z)).

We define ψ = (g -i)/(g + i).

Now, for n ≥ 2 and 1/16π(n + 1) ≤ h < 1/16πn, one gets hence:

= ω Ω (a, {w ∈ ∂Ω ; Im w > 4πn}) ≤ ω Ω (a, ∂Ω \ ∂Ω n ) ≤ δ n ≤ δ(h) , proving Lemma 5.5.


A Gramfort

G Varoquaux

B Thirion

Beyond brain reading: randomized sparsity and clustering to simultaneously predict and identify

The prediction of behavioral covariates from functional MRI (fMRI) is known as brain reading. From a statistical standpoint, this challenge is a supervised learning task. The ability to predict cognitive states from new data gives a model selection criterion: prediction accuracy. While a good prediction score implies that some of the voxels used by the classifier are relevant, one cannot state that these voxels form the brain regions involved in the cognitive task. The best predictive model may have selected by chance non-informative regions, and neglected relevant regions that provide duplicate information. In this contribution, we address the support identification problem. The proposed approach relies on randomization techniques which have been proved to be consistent for support recovery. To account for the spatial correlations between voxels, our approach makes use of a spatially constrained hierarchical clustering algorithm. Results are provided on simulations and a visual experiment.

Introduction

Functional MRI (fMRI) is an imaging technique that measures Blood Oxygen-Level Dependent signal changes caused by brain activity. Detecting and localizing these changes can be used to improve our understanding of brain function. Over the last decade, many contributions have proposed to tackle this challenge using statistical learning and more specifically supervised learning methods [START_REF] Mitchell | Learning to decode cognitive states from brain images[END_REF]. The data are fMRI volumes -3D images made of voxels-and the target to predict is, for example, the stimulus that was presented to the subject in the scanner. This formulation of the problem is commonly called brain reading or decoding.

A strong benefit that supervised learning methods bring to brain mapping is the ability of the estimator to account for a distributed pattern of active voxels. While standard statistics for brain mapping model only one voxel at a time, or local clusters, brain reading can be applied to full brain data. The method is said to be multivariate: the learned prediction function relies on correlations between distant brain regions.

To actually achieve brain mapping, the learning method used for decoding should inform about which voxels are useful for the prediction. This constraint naturally favors linear classifiers for which the prediction function is obtained from a linear combination of the voxel amplitudes. We call the coefficients of this linear combination the weights of the estimator. They form a spatial map.

Functional MRI data are a spatially smoothed representation of the underlying neural signals. Consequently, the activations are not only distributed over the entire brain but also spatially correlated. For better prediction performance, the estimators should incorporate this prior knowledge. It is natural to promote prediction functions relying on only a few brain regions, for instance using sparsity inducing regularization methods [START_REF] Krishnapuram | Sparse multinomial logistic regression: Fast algorithms and generalization bounds[END_REF][START_REF] Yamashita | Sparse estimation automatically selects voxels relevant for the decoding of fMRI activity patterns[END_REF][START_REF] Carroll | Prediction and interpretation of distributed neural activity with sparse models[END_REF]. In addition, to account for the spatial structure in the signal, the estimator can make use of the three-dimensional grid structure over which the signal is defined. This can be achieved with convex regularization promoting piecewise constant weights [START_REF] Michel | Total variation regularization for fMRI-based prediction of behavior[END_REF] or by constructing hierarchically organized spatial features using a spatially constrained Ward clustering method [START_REF] Ward | Hierarchical grouping to optimize an objective function[END_REF] and learning a linear decision function defined over this new set of features. [START_REF] Michel | A supervised clustering approach for fMRI-based inference of brain states[END_REF] perform the learning step with a greedy top-down approach while [START_REF] Jenatton | Multi-scale Mining of fMRI data with Hierarchical Structured Sparsity[END_REF] use a hierarchical convex prior.

A caveat with the decoding approach is that the model is selected to optimize the prediction, while the localization of brain function requires instead to optimize the identification of the brain regions involved in the task. A good prediction indicates that the identified regions are sufficient to predict but it means neither that they are the true ones nor that they cannot be better estimated. A simple illustration is that different weights can lead to the same prediction accuracy [START_REF] Rish | Sparse regression analysis of task-relevant information distribution in the brain[END_REF]. In this paper we address the identification problem while taking into account the specificity of fMRI data: distributed patterns and spatial correlations.

Notations We write vectors with bold letters, a ∈ R n , matrices with capital bold letters, A ∈ R n×n . We denote

a 1 = n i=1 |a i | the 1 norm and a 2 = n i=1 a 2 i the 2 norm.
2 Method for sparse recovery with spatial structure

Let us consider the linear classification model:

y = sign(x w + b) = sign p v=1 x v w v + b , (1)
where y ∈ {-1, 1} represents the target to predict, sign stands for the sign function, p is the number of voxels in the grid, and (w, b) are the model parameters to be estimated: the weight vector and the intercept, also called bias term. The vector x ∈ R p is an fMRI activation volume. Using a logistic regression model, the estimation problem over a training set formed by n volumes reads:

( w, b) = argmin (w,b) 1 n n i=1 log 1 + e -y i (x i w+b) + λΩ(w), λ > 0, (2)
where λ controls the level of regularization and Ω is the regularizing function that typically promotes sparse and potentially spatially structured weights. Sparse logistic regression (SLR) refers to the case where Ω(w) = w 1 . We denote by X ∈ R n×p the design matrix formed by the concatenation of all x i .

One of the issues with sparse methods is the instability of the estimated support of w, particularly when the columns of X are very correlated as it is the case with fMRI data. To stabilize the estimated support, it has been proposed to randomly perturb the design matrix [START_REF] Meinshausen | Stability selection[END_REF] by taking only a fraction of the training samples and randomly scale each column, in our case each voxel. By repeating the later procedure and then counting how often each voxel is selected across the repetitions, each voxel can be assigned a score. The higher is the score, the more likely is the voxel likely to be truly informative. In a regression setting, this procedure is called Randomized Lasso [START_REF] Meinshausen | Stability selection[END_REF]. We apply it here to a binary classification problem (see e.g. [START_REF] Ng | Feature selection, l1 vs. l2 regularization, and rotational invariance[END_REF]).

Let k ∈ [1 . . . K] denote the repetition and w k be the corresponding estimated weight vector. The design matrix X k is formed by a random fraction π of the training data. Each column of X k is then randomly scaled to 1 or to 1 -a with equal probability. The procedure is a subsampling of the data and a random perturbation of each voxel. The stability score of each voxel v is then the percentage of the repetitions for which the voxel has a non-zero weight, i.e., is used for the prediction. A voxel v is used if the corresponding entry in the weight vector w k estimated at repetition k is non-zero. We denote it by v ∈ supp(w k ). The stability score can then defined as

s v = #{k s.t. v ∈ supp(w k )}/K ∈ [0, 1]
. The estimated support is then defined as {v s.t. s v ≥ τ }. In the following experiments a is set to 0.5 and π to 75%, while τ is estimated by cross-validation in a discrete set of values T . Following experiments use a fixed value T = {0.25} or a grid T = {0.1, 0.2, 0.3, 0.4, 0.5}. For every threshold τ ∈ T , a cross-validation score is obtained using a 2 -logistic regression model (Ω(w) = w 2 2 ). The estimated τ is the one yielding the highest score.

To improve the stability of the estimation and inform the estimator about the grid structure of the data, we propose using Ward hierarchical clustering as in [START_REF] Michel | A supervised clustering approach for fMRI-based inference of brain states[END_REF] to train the classifiers on data-driven spatial clusters. At each repetition the brain volume is first segmented in q spatially connected regions in which the fMRI signal is averaged. As the hierarchical tree is estimated each time on a random fraction of the data, the tree is different for every randomization. Note that a similar procedure is performed in the Random Forests algorithm [START_REF] Breiman | Random forests[END_REF]. One obvious benefit of this procedure is that it tends to produce an "average" tree which balances with the greedy hierarchical construction of a single tree. The SLR is then fitted on a q-dimensional dataset. A voxel is marked as active in repetition k if it belongs to a region with a non-zero weight. Although the estimated supp(w k ) is in R q , we will still write v ∈ supp(w k ).

The main benefit of the additional clustering step is to reduce the correlations in the design matrix, therefore improving the behavior of sparse methods. Our method can thus select more voxels than the number of observations, which would be impossible with standard SLR and difficult with only randomization.

The procedure is summarized in Algorithm 1.

Algorithm 1 Randomized Sparse Logistic Regression with hierarchical features

Input: Set 0 < a < 1, T (e.g. {0.1, 0.2, 0.3, 0.4, 0.5}), K (e.g. 200), π (e.g. 0.75).

1: Estimate q and λ with cross-validation 2: for k = 1, ..., K do 3:

Randomize design X k with data subsampling and random feature scaling 4:

Hierarchical clustering (segment brain in q regions) 5:

Estimate

w k ∈ R q with SLR (

Simulation study and fMRI results

We now present a simulation study followed by results on fMRI data recorded during an object recognition experiment. Experiments were performed with the scikit-learn [START_REF] Pedregosa | Scikit-learn: Machine learning in Python[END_REF] using LibLinear [START_REF] Fan | LIBLINEAR: A library for large linear classification[END_REF].

Simulation study -The simulation data consist of training and test sets each of 160 volumes. Each volume is a cube of (9×9×9 voxels). The active regions are 2 cubes of size 2×2×2 located at two opposite sides of the cube (see Fig. 1-a). Data are corrupted by a Gaussian additive noise and smoothed. The parameters λ and q are estimated by 5-fold cross-validation on the training set. Then stability scores are estimated with K = 200 repetitions. Figure 1 presents the F-values for each voxel, as in conventional brain mapping, and the selection scores s. Accuracy is quantified, for the identification, with a Precision-Recall (PR) curve on the recovered support and, for the prediction, with a Receiver-Operating-characteristic (ROC) curve on the predicted labels. Prediction performances using the known true weights are also given as baseline. A first interesting observation is that although SLR outperforms a linear SVM and a 2 -logistic regression for prediction (Figure 1-e), it is clearly worse for the identification (Figure 1-d). This illustrates that the model that predicts the best may not be the model built from the true active voxels. What is also interesting it that the proposed method clearly outperforms all alternative methods for support recovery, while also giving almost optimal prediction accuracy. Identification results with different active regions sizes are present in Tab. 1. Our approach consistently provides the best estimation, except when the solution is very sparse (1 voxel), in which case it is outperformed by sparse estimators. fMRI data analysis -The investigated fMRI data consist of five subjects recorded during a visual object recognition protocol [START_REF] Haxby | Distributed and overlapping representations of faces and objects in ventral temporal cortex[END_REF]. In this experiment subjects were asked to recognize 8 different types of objects. We focus here on the binary classification task that consist in predicting whether the subject was viewing a house or a face. The data consist of 12 sessions, that were split in a training set and an independent test set. Each session contained 18 volumes (9 in each category). Preprocessing of the data consisted in motion correction using FSL MCFLIRT and a removal of linear trends in each session. We first present results using respectively 4 and 8 sessions for the train and test sets. The first subject is presented in details in Fig. 2. The ROC curves show that the best prediction accuracy is obtained with the proposed method, followed by 1 -logistic regression and then the 2 penalized methods (Logistic and linear SVM). As shown on Fig. 2-a, voxels with strong selection scores (above 0.25) are located within Fusiform gyrus in a region known as the Fusiform Face Area (FFA). The ROC curves for the other subjects are presented in Fig. 3. The mean ROC area across subjects is 0.989 while it is only 0.869 for the SLR, 0.807 for the 2 -LR and 0.808 for the linear SVM. These results show that the proposed method consistently outperforms alternative approaches in terms of prediction accuracy. This method also yields a spatially structured and a meaningful estimated support in the FFA which suggests that the randomization procedure employed improves the support recovery as shown in [START_REF] Meinshausen | Stability selection[END_REF].

In order to further investigate the performance of the method as a function of the number of training data, we have conducted the same experiments when varying the number of sessions used from estimating the support and fitting the predictive model. Results are presented in Figure 4. A first interesting observation is that all methods tend to predict almost perfectly when using a large training set. Another observation is that here again, a linear SVM and an 2 -logistic regression yield very similar results. The SLR outperforms the later methods when using more than 3 sessions for fitting the model. Finally, the proposed method is the only one yielding almost predictions as soon as the number of sessions exceeds three.

Conclusion

In this work, we have shown that a randomization technique coupled with a spatial clustering algorithm could significantly improve the identification of predictive brain regions while yielding better prediction scores. The sparse randomized logistic used for that purpose allowed to stabilize the support estimation while the clustering pre-processing addresses successfully the problem of strong spatial correlations. This contribution illustrates a somehow unintuitive fact that among the set of models, like the one obtained with a sparse method when varying the regularization parameter, the model that predicts the best is not always the model that identifies best the good voxels. The optimization of the prediction score on unseen data or the identification of the good voxels can lead to different models. The nice observations presented in this work, is that the proposed procedure improves both aspects, the support identification and the prediction scores.

Fig. 1 .

 1 Fig. 1. Simulation results: a. ground truth, b. F-values, c. scores, d. identification precision-recall, e. prediction ROC

Fig. 2 .

 2 Fig. 2. Results on fMRI object recognition task (face vs. house) on subject 1. a. Scores with Ward Randomized Logistic Regression on the selected voxels in the Fusiform Face b. The red isosurface shows the voxels with scores above 0.25. The green surface delimits the brain volume in which the fMRI data are defined. c. Prediction ROC.

Fig. 4 .

 4 Fig. 4. ROC AUC scores as a function of the number of training sesssions.

Table 1 .

 1 Area

		2)
	6: end for
	7: Set scores sv = #{k s.t. v ∈ supp(w k )}/K ∈ [0, 1]
	8: Set estimated support {v s.t. sv ≥ τ } (τ ∈ T estimated by cross-validation with
	2-logistic regression)
	under	Methods Ward Rand. LR F-test 1-LR 2-LR SVM
	Precision-Recall curve	1×1×1	0.84	0.589 1.0 0.773 0.773
	as a function of the	2×2×2	0.98	0.829 0.523 0.775 0.775
	active region size.	3×3×3	0.786	0.749 0.456 0.535 0.631


Keywords:

Localisation de sources sonores en utilisant l'estimation de delais temporels contraints gometriquement.

Résumé : Dans cet article nous prsentons une mthode de temps gomtriquement contraints retard d'estimation pour la source sonore la localisation (gTDE). Une analyse algbrique rvle que la mthode peut traiter avec un nombre arbitraire de non-coplanaires microphones. Nous tirons un problme d'optimisation sous contrainte non-linaire qui peut łtre rsolu en utilisant la programmation convexe. contrairement techniques existantes, qui considrent paires TDE, la mthode propose estime de manire optimale un ensemble de retards qui sont compatibles avec l'emplacement de la source. De vastes expriences simules de valider la mthode en prsence de bruit et de rverbrations.

Introduction

For the last two decades, source localization from time delay estimates (TDE's) has proven to be an extremely useful methodology with a variety of applications in such diverse fields as aeronautics, telecommunications and robotics. Also referred to as multilateration, this problem is highly related to the one of estimating time delays. We are particularly interested in the development of a general-purpose TDE-based method for sound-source localization in indoor environments, e.g, human-robot interaction, ad-hoc teleconferencing microphone arrays, etc. This type of consumer-oriented applications are extremely challenging for several reasons: (i) there may be several sound sources and their number varies over time, (ii) regular rooms are echoic, thus leading to reverberations, and (iii) the microphones are often embedded in devices (robot heads, smart phones, etc.) generating high-level noise.

Structure 2 Related work

The TDE problem was very well investigated: A recent review can be found in [START_REF] Chen | Time Delay Estimation in Room Acoustic Environments: An Overview[END_REF]. The vast majority of existing approaches deals with a microphone pair but it is not straightforward to extend most of these methods to several microphones. Methods addressing multichannel TDE can be roughly divided into two categories: methods estimating the acoustic impulse responses and methods exploiting the redundancy among several microphones. [START_REF] Doclo | Robust Adaptive Time Delay Estimation for Speaker Localization in Noisy and Reverberant Acoustic Environments[END_REF] is illustrative of the first category: a method based on generalized eigenvalue decomposition is proposed. The second category is represented by [START_REF] Chen | Robust time delay estimation exploiting redundancy among multiple microphones[END_REF] where a multichannel criterion based on cross-correlation is proposed to estimate time delays using a linear microphone array. In both cases, experiments are performed on speech data in a simulated indoor environment.

As already mentioned, an alternative to TDE is multilateration, which makes assumptions about the time delay estimates. This provides a framework for casting the problem into maximum-likelihood estimation or into mean-squared error minimization (see [START_REF] Pertilä | Acoustic Source Localization in a Room Environment and at Moderate Distances[END_REF] for a review). Two recent methods deserve to be mentioned. In [START_REF] Brutti | Comparison between different sound source localization techniques[END_REF] the authors use the acoustic maps together with the GCC-PHAT technique to localize sound sources from TDE's. The model in [START_REF] Ribeiro | Using reverberation to improve range and elevation discrimination for small array sound source localization. Audio, Speech, and Language Processing[END_REF] includes the reverberations in order to enhance the localization performance while using a uniform circular array of microphones.

There are, of course, tones of work that has been done before, related to the topic. I will briefly describe several works. In the next section, some selected papers will be described in more detail.

• The method described in [START_REF] Chan | A simple and efficient estimator for hyperbolic location[END_REF] performs 2D localization from a non-linear array of three microphones. The use an unbiased estimator for TDOA. From this paper, I derived the 3D localization from a non-coplanar array of four microphones. This is not a contribution, since the analogy with the 2D is straightforward and too obvious.

INRIA

• In [START_REF] Smith | Closed-form least-squares source location estimation from range-difference measurements[END_REF][START_REF] Friedlander | A passive localization algorithm and its accuracy analysis[END_REF] they derive a least-squares method for more than four microphones assuming independent noise in the TDOA measurements.

• In [START_REF] Brandstein | A closed-form location estimator for use with room environment microphone arrays[END_REF][START_REF] Brandstein | A practical methodology for speech source localization with microphone arrays[END_REF] they approximate the hyperboloids by cones and perform linear intersection in the least-squares sense. Tested on simulated and real data respectively. TDOAs are estimated independently.

• A linear correction for least-squares approach is presented in [START_REF] Huang | Real-time passive source localization: A practical linear-correction least-squares approach[END_REF]. Just on simulated data.

• An outlier detection algorithm from TDOA noisi measurements by parity check is presented in [START_REF] Galati | New Approaches to Multilateration processing: analysis and field evaluation[END_REF].

• Numerical relaxation methods for the hard-least-squares problem are presented in [START_REF] Beck | Exact and approximate solutions of source localization problems[END_REF] to provide an approximate solution. Just on simulations.

• A method is presented for non-line-of-sight TDOA localization in [START_REF] Chan | Time-of-arrival based localization under NLOS conditions[END_REF].

• A method is presented based on TOA (not TDOA) in [START_REF] Cheung | Least squares algorithms for timeof-arrival-based mobile location[END_REF].

• Two methods are presented in [START_REF] Doclo | Robust Adaptive Time Delay Estimation for Speaker Localization in Noisy and Reverberant Acoustic Environments[END_REF] to extend the adaptive eigenvalue decomposition (see [START_REF] Chen | Time Delay Estimation in Room Acoustic Environments: An Overview[END_REF] for a nice review on methods): stochastic gradient and data prewhitening. Simulated data.

• In [START_REF] Drake | Geolocation by time difference of arrival using hyperbolic asymptotes[END_REF] the problem of localization from TDOA with more than four microphones is casted into a first order Taylor series approximation, which geometrically corresponds to a hyperbolic asymptote localiation algorithm. Simulated data.

• A really nice survey on time difference estimation (TDE = TDOA) methods is presented in [START_REF] Chen | Time Delay Estimation in Room Acoustic Environments: An Overview[END_REF]. The comparisons are run extensively. Two microphone, Nmicrophone and adaptive methods are presented.

• A method for multi-microphone TDOA estimation is presented in [START_REF] Chen | Robust time delay estimation exploiting redundancy among multiple microphones[END_REF]. Simulated data. Linear array.

• Here a multichannel cross-correlation coefficient os presented (applied to TDE with a linear array). Varechoich chamber.

• A method for ML estimation in the near-field and cross-bearing estimation from separated arrays is presented in [START_REF] Chen | Acoustic source localization and beamforming: theory and practice[END_REF]. Single/Multiple source. Real experiments in a semi-anechouc chamber. 2D localization.

• In [START_REF] Handzel | Biomimetic sound-source localization[END_REF][START_REF] Handzel | High Acuity Sound-Source Localization by means of a Triangular Spherical Array[END_REF] a study in physics on how waves are refracted.

• Localization from ILD in [START_REF] Hinich | Localization of Sound Source in 3D Space[END_REF].

• TDOA & FDOA based localization and velocity estimation of sound sources in [START_REF] Ho | Source localization using TDOA and FDOA measurements in the presence of receiver location errors: Analysis and solution[END_REF].

• Humanoid binaural sound trakcing using kalman filtering and HRTF's (simulated and real data) in [START_REF] Keyrouz | An Enhanced Binaural 3D Sound Localization Algorithm[END_REF][START_REF] Keyrouz | Humanoid Binaural Sound Tracking Using Kalman Filtering and HRTFs[END_REF].

• A very interesting paper on kalman filters for speaker tracking using real data in [START_REF] Klee | Kalman Filters for Time Delay of Arrival-Based Source Localization[END_REF].

• There is a theoretical study of the close-dorm solution to a really limited (almost ridiculous) extent in [START_REF] Mellen | Closed-form solution for determining emitter location using time difference of arrival measurements[END_REF].

• A method for sound source localization using ILD and ITD with less than four microphones is presented in [START_REF] Sam Ge | Robust Sound Localization using Small Number of Microphones[END_REF].

• A nice survey of mathematical methods for indoor localization is presented in [START_REF] Seco | A survey of mathematical methods for indoor localization[END_REF].

• ML framework for localization in [START_REF] Sheng | Maximum likelihood multiple-source localization using acoustic energy measurements with wireless sensor networks[END_REF]. Distributed sensors. TDOA errors are independent and Gaussian.

• CRLB optimal mapping from all TDOAs to N-1 TDOAs. TDOA errors are independent and Gaussian. Simulated data. ( [START_REF] So | Closed-form formulae for time-difference-of-arrival estimation[END_REF])

• [START_REF] Strobel | Classification of time delay estimates for robust speaker localization[END_REF] proposes a ML method to classify estimated TDOAs based on learning conditional pdfs. Tested with real data from an anechoic chamber.

• In [START_REF] Urruela | Novel closed-form ML position estimator for hyperbolic location[END_REF] the authors propose a new ML estimator (I did not get what is new on it). The TDOA measurements are independent. Simulated data.

• Optimal microphone arrays design to minimize the CRB (interesting). [START_REF] Yang | Cramer-Rao bound and optimum sensor array for source localization from time differences of arrival[END_REF].

• Theoretical study of 2D sensor arrays in [START_REF] Yang | A theoretical analysis of 2D sensor arrays for TDOA based localization[END_REF].

• Relaxation methods for localization using TDOAs and simulated data in [START_REF] Yang | Efficient convex relaxation methods for robust target localization by a sensor network using time differences of arrivals[END_REF].

• A ML framework for localization using directional microphones is presented in [START_REF] Zhang | Maximum Likelihood Sound Source Localization and Beamforming for Directional Microphone Arrays in Distributed Meetings[END_REF][START_REF] Zhang | Maximum Likelihood Sound Source Localization for Multiple Directional Microphones[END_REF] and tested with real data.

• Modification of the PHAT algorithm to discard frequencies in which the signal is not present. Localization and tracking with real data in [START_REF] Valin | Robust 3D localization and tracking of sound sources using beamforming and particle filtering[END_REF].

Closest work

There are some papers among the previous list which are more related to what I would like to do:

1. The method presented in [START_REF] Chen | Robust time delay estimation exploiting redundancy among multiple microphones[END_REF] derives a multi-microphone TDOA estimator. That is, it estimates all the TDOA values at ones, instead of each TDOA separately. It has been just tested on simulated data!!!! 2. The method in [START_REF] Chen | Acoustic source localization and beamforming: theory and practice[END_REF] presents a method for sigle/multiple source TDOA estimator. This is interesting, because it is not just the single source case.

3. The method in [START_REF] Valin | Robust 3D localization and tracking of sound sources using beamforming and particle filtering[END_REF] modifies the PHAT algorithm to discard frequencies in which the signal is not present.

INRIA 3 Signal model

In this section, the signal acquisition model is described. The most important elements of this model are: the sound source position (S), the number of microphones (M ) and their positions (M 1 , . . . , M M ) and the emited signal (x(t)).

From this basic elements, we can write the expression for the received signal at the m-th microphone:

x m (t) = x(t -t m ) + n m (t), (1)
where n m is the microphone's noise and t m is the time-of-arrival. The physical model for t m assumes constant sound propagation speed, denoted by ν. Hence we write t m = S -M m /ν. Using this model, the expression for the time dealy between the m-th and the n-th microphones arises naturally. Denoted by t m,n , this TDE is generated by the following geometric model:

t m,n = t n -t m = S -M n -S -M m ν . (2)
Notice that, for a fixed value of t m,n and given the microphones' positions, the sound source generating t m,n lies in a two-sheet hyperboloid with focci M n and M m . It is very important to notice also, that the set of all M (M -1)/2 time delay values corresponding to the M microphones are not independent, since t m,n = t m,k + t k,n .

In the following section we develop a method for time delay estimation from the signal model described in (1).

Multichannel Time Delay Estimation

In this section we present the framework used to estimate the time delay between microphone signals. As in [START_REF] Chen | Robust time delay estimation exploiting redundancy among multiple microphones[END_REF], the signal acquired by one of the microphones is set as the reference signal. The other signals are properly scaled and delayed in order to predict the reference signal. The best prediction, that is the best scaling and delaying values, correpond to the time delay estimates we are looking for. More precisely, given the signals acquired by the M microphones, x 1 (t), . . . , x M (t), we can set x 1 as the reference signal and write the following linear prediction error:

e c,t (t) = x 1 (t) - M m=2 c 1,m x m (t + t 1,m ), (3)
where c = (c 1,2 , . . . , c

(c * , t * ) = arg min c,t J(c, t) = arg min c,t E e 2 c,t (t) . (4)
If we develop the expression above and derive with respect to c, we can show that c * (t) = R -1 (t)r(t), with:

R(t) =     R 2,2 (0) R 2,3 (t 1,3 -t 1,2 ) • • • R 2,3 (t 1,3 -t 1,2 ) R 3,3 (0) • • • . . . . . . . . .     , r(t) =     R 1,2 (-t 1,2 ) R 1,3 (-t 1,3 ) . . .     . In the expressions above R i,j (τ ) = E {x i (t)x j (t -τ )}. Notice that R i,j (τ )
has its maxima at the true value of t i,j . In some sense we are looking for the combination of t 1,k that maximize the criteria, that turns out to maximize the cross-correlation functions. The optimization problem is:

t * = arg min t J(c * (t), t) = arg min t R 1,1 (0) -r T R -1 r , (5)
where we suppressed the explicit dependency of r and R on t. Also, it can be shown (see Appendix A) that this problem is equivalent to the following one:

t * = arg min t J(t) = arg min t det R , (6)
where R is the matrix of normalized cross-correlation functions evaluated at t:

R =      1 ρ 1,2 (-t 1,2 ) ρ 1,3 (-t 1,3 ) • • • ρ 1,2 (-t 1,2 ) 1 ρ 2,3 (t 1,3 -t 1,2 ) • • • ρ 1,3 (-t 1,3 ) ρ 2,3 (t 1,3 -t 1,2 ) 1 • • • . . . . . . . . . . . .      . ( 7
)
As in [START_REF] Chen | Robust time delay estimation exploiting redundancy among multiple microphones[END_REF] the time delay estimation problem is cast into a non-linear optimisation problem. However, since no assumption has been done about the geometry of the microphone array, the problem is multidimenional (M -1 variables) instead of unidimensional.

Geometric model

Depending on the geometry of the microphone array, some values of t do not correspond to a position in the sound source space. This was fully exploited by Chen in the particular case of linear arrays by constraining the signal model. Since this is not possible in general, we need a change of paradigm when dealing with arbitraty geometry arrays. Hence, we propose to modify the estimation procedure by adding a non-linear constraint. This constraint describes the geometric feasibility of the values in t. That is why the proposed algorithm is called geometrically-constrained time delay estimation (gTDE).

In this section we present the geometric model for TDE, given the sound source's position. We will call that model the "Direct". After we derive a method to inverse this direct model, i.e., to obtain the sound source's position from the TDEs. We call that the "Inverse" model.

INRIA

Direct model

Given the microphones' position, M 1 , . . . , M M , and the sound source's position, S, the geometric generative model for the TDEs is:

t ij (S) = S -M j -S -M i ν , ( 8
)
where ν is the sound speed. We can also think about the following application:

F : R 3 → R M (M -1)/2 S → F (S) := (t ij (S)) i<j . (9)
Notice that the M (M -1)/2 TDEs are not independent. Actually the following relation holds

t ij = t kj -t ki . (10)
From this relation we observe that the maximal subsets of independent TDEs are of size M -1. In other words, all subsets of size M have at least one TDE which depends on the other M -1. For example, we can think on the subset of all TDEs related to one microphone: {t k1 , . . . , t kk-1 , t kk+1 , . . . , t kM }.

Inverse model

The model presented before is invertible under some hypothesis. This means that we can recover the 3D position of the sound source (S) from a set of TDEs. In the next we will derive this inverse geometric model. As explained in the previous section we can use just M -1 TDEs since the other are dependent on them. We will use the TDEs related tot he first microphone: {t 12 , . . . , t 1M }. To be precise we will use a scaled version of the TDE for the geometric derivation: the DDE (Distance Difference of Arrival), noted by {d 12 , . . . , d 1M }, where

d ij = ν t ij .
Let us denote by d m (S) the distance from the sound source S to the m-th microphone, which can be written as:

d 2 m (S) = S -M m , S -M m = S 2 -2 S, M m + M m 2 . ( 11
)
Using the fact that d m (S) = d 1m (S) + d 1 (S), we can write:

d 2 1m (S) + 2d 1m (S)d 1 (S) + d 2 1 (S) = S 2 -2 S, M m + M m 2 , (12)
and subtract the following equation:

d 2 1 (S) = S 2 -2 S, M 1 + M 1 2 . ( 13
)
This leads to:

d 2 1m + 2d 1m d 1 = M m 2 -M 1 2 -2 S, M m -M 1 , (14)
where we suppressed the explicit dependency of d m1 and d 1 on S. Since this equation is non-trivial for m ≥ 1, we can write the following system of equations:

           d 2 12 + 2d 12 d 1 = M 2 2 -M 1 2 -2 S, M 2 -M 1 . . . d 2 1M + 2d 1M d 1 = M M 2 -M 1 2 -2 S, M M -M 1 (15)
Notice that the system of equations is linear in S. Indeed, we can rewrite it in the following form:

M • S = K -D 2 -2D d 1 (16)
where

M = 2     (M 2 -M 1 ) T . . . (M M -M 1 ) T     , K =    M 2 2 -M 1 2 . . . M M 2 -M 1 2    , D =    d 12 . . . d 1M    , (17)
and D 2 is the component-wise square power of D.

Notice that the matrix M is full rank if and only if the M microphones do not lie all in the same plane. This means that the problem has non-degenerated solution when the microphones are non co-planar. We will refer to that as the non-degeneracy condition of the problem. Assuming that the non-degeneracy condition is satisfied, we can solve the system of equations (linear in S), by inversion:

S = M ‡ K -D 2 -2M ‡ D d 1 = Ad 1 + B, (18)
where A = -2M ‡ D and B = M ‡ K -D 2 , and M ‡ is the pseudoinverse matrix of M. We now substitute this solution inside the equation for d 1 :

d 2 1 = S -M 1 , S -M 1 = Ad 1 + B -M 1 , Ad 1 + B -M 1 = A 2 d 2 1 + 2 A, B -M 1 d 1 + B -M 1 2 . ( 19
)
This leads to a second order equation for d 1 which we can solve analytically and substitute in the equation for S.

d + 1 = A, B -M 1 + √ ∆ A 2 -1 d - 1 = A, B -M 1 - √ ∆ A 2 -1 , (20)
where

∆ = A, B -M 1 2 -B -M 1 2 ( A 2 -1). (21)
The two solutions are:

S + = Ad + 1 + B S + = Ad - 1 + B. (22)
INRIA Since D = νt, we can rewrite ∆ depending on t: ∆(t). Notice that the condition for the existence of a solution is ∆ ≥ 0; in such case we will call the set of TDE consistent. Also if ∆ > 0 two solutions exist. However, just one of them is correct. Actually, we observed that the following statement holds:

F (S + ) = F (S) F (S -) = -F (S) or F (S + ) = -F (S) F (S -) = F (S)
Geometrically speaking, we are looking fot the intersection of M -1 sheets of M -1 two-sheet hyperboloids. These intersect just in one point. Notice, however, that we can change the sign of the TDEs and we will get the same solutions to the equation. This means that we are considering also the intersection of the other M -1 one-sheet hyperboloids (the other sheets of the two-sheet hyperboloids). This new intersection provides us the second and fake solution.

Geometrically constrained criterion

We would like, hence, to put together the signal prediction criterion subject to the geometric criterion. The optimization problem is set as:

(t * ) = arg min t J(t) s.t. ∆ (t) ≥ 0. ( 23
)

Implementation, Experiments & results

The minimization of ( 23) is carried out using a publicly available MATLAB implementation [START_REF] Carbonetto | MATLAB primal-dual interior-point solver for convex programs with constraints[END_REF] of the log-barrier interior point method [START_REF] Boyd | Convex Optimization[END_REF]. This method is designed for continuous convex optimization problems. On one hand, it is likely to fail in finding the global optimum of non-convex problems such as [START_REF] Keyrouz | An Enhanced Binaural 3D Sound Localization Algorithm[END_REF]. To overcome this issue, the algorithm starts from several initial points, i.e., the set S I = {t I i } P i=1 . For each one of these initializations a local minimum is found, then the minimum over these local minima is selected. In our simulations, P = 4096 and S I consists of points placed in a regular rectangular grid. Since the t m,n 's have upper and lower bounds ( M m -M n and -M m -M n respectively), the grid limits are defined by the geometry of the problem. On the other hand, the function to optimize must be continuous and the signals are discrete. We chose to compute the normalized cross-correlation function of the linear interpolation of the discrete signals.

In order to accurately evaluate and validate the proposed gTDE method, we developed a formal evaluation protocol using simulated data. A 3×4×2.5 meter room, with uniform absorption coefficients, was simulated using the state-of-the-art Image-Source Model (ISM), [START_REF] Lehmann | Prediction of energy decay in room impulse responses simulated with an image-source model[END_REF] available at [START_REF] Lehmann | Matlab code for image-source model in room acoustics[END_REF]. The main parameter of this model is T 60 , which The source emitted speech fragments randomly chosen from [START_REF] Garofolo | Timit acoustic-phonetic continuous speech corpus[END_REF]. One hundred millisecond cuts of these sounds are the input of the evaluated method. We assume that only one source is emitting during these cuts. Finally the sensor noise, whose power depends on the chosen SNR, is added to these cuts.

Fig. 1 shows the results obtained with four different methods: pair-wise independent estimation of the time delays (bypairs), estimation using multichannel information without minimization (init), estimation based on unconstrained time delay minimization of (6) (tde), and the proposed geometrically-constrained minimization (gtde). Fig. 1(a) plots the percentage of good (non-anomalous) estimates, i.e., with an absolute error smaller than 100 µs and Fig. 1(b) plots the standard deviation of the good estimates as a function of the signal-to-noise ratio. These results correspond to an anechoic setup, T 60 = 0. As expected, the proposed method significantly improves the percentage of good TDEs while it lowers down the standard deviation.

Additional simulations were carried out in order to precisely evaluate the performance of the gTDE method in the presence of noise and of reverberations. Figures 2 and3 show the results on time delay estimation and sound source localization for different levels of noise and reverberations. In all the plots, the x-axis represents the SNR value (dB). The color corresponds to the method used: green for bypairs and blue for gtde; and the line style corresponds to the level of reverberation: solid-circle for T 60 = 0 s and dashed-square for T 60 = 0.1 s.

Regarding the TDE, in Fig. 2 clearly outperforms the baseline. This is not a surprise since it is using all the available information to consistently estimate the TDEs at a time. A remarkable fact is that the proposed method under reverberant conditions has similar performance than the baseline method in the anechoic case. Also, with higher percentage of non-anomalous estimates, the gtde method has lower error standard deviation.

Concerning the localization, Fig. 3(a) plots the percentage of localization inliers and Fig. 3(b) the standard deviation of the angular error. A sound source is considered to be an inlier if the angular absolute error is less than 30 • . As in the case of the TDE, the methods' performance improve with the SNR. Notice also that the proposed method outperforms the baseline, even when the proposed method is under echoic conditions and the baseline is under anechoic conditions. Generally speaking, the methods perform as expected. The higher the SNR value the better the methods estimate the time delays, the higher the percentage of inliers and the lower the localization error. We can also see a clear trend with respect to the reverberation level: the methods' performance decreases with T 60 . However this two parameters have different effects on the function to minimize. On one side, the sensor noise decorrelates the microphones' signals leading to much more (and randomly spread) local minima and increasing the value of the true minimum. If this effect is extreme, the hope for a good estimate decreases fast. On the other side, the reverberations produce a little amount of strong local minima. This perturbation is systematic given the source position in the room. Hence, there is hope to learn the effect of such reverberations in order to improve the quality of the estimates. These two types of perturbations of the function to minimize clearly have different effects on the results. Notice that the reverberation level has almost no effect on the quality of the estimates when the SNR is low. However, when this random effect disappears (i.e., higher values of SNR) a systematic and significative difference appears both in the time delay and in the localization estimates.

Finally, the authors would like to remark on the method's robustness. Notice that, for moderate levels of noise and reverberations (T 60 ≤ 100 ms, SNR ≥ 0 dB) the method is able to localize the sound source with mean angular squared error of 6 • in more than 60 % of the cases.

Conclusions & future work

In this paper, a new method on time delay estimation for sound source localization working on non-coplanar microphone arrays is presented. The estimation is cast into a multivariate optimization problem. In addition, the geometric model for the time delays add a non-linear constraint. The optimal value of the problem is a consistent set of time delays, useful to localize the sound source. Experiments on simulated data show the quality of the method and validate the approach.

The are several ways to extend this work. As outlined before, it would be very useful to learn the effect the reverberations have on the objective function as in [START_REF] Ribeiro | Using reverberation to improve range and elevation discrimination for small array sound source localization. Audio, Speech, and Language Processing[END_REF]. Also, it is worth to consider the multiple source case, following approaches like [START_REF] Chen | Acoustic source localization and beamforming: theory and practice[END_REF]. Besides that, a frequency decomposition stage may be useful to avoid the analysis in non-informative frequency bands ( [START_REF] Valin | Robust 3D localization and tracking of sound sources using beamforming and particle filtering[END_REF]). Also, experiments on more reverberant data, and on real data have to be done in order to explore the real extent of these initial results. Last but not least, it would be interesting to explore cases in which the microphones' positions have some error, and see how to adapt the method to estimate the time delays and correct the microphones' positions.

A Criteria equivalence

This appendix is devoted to prove the equivalence of both optimization problems: (5) and ( 6). We will start from the expression of the second criterion and get to the expres-INRIA sion of the first (original) one. The objective function to minimize is:

det R = det (R) - M i=2 (-1) i+1 ρ 1,i det (R i ) ,
where

R i =   r R 1 • • • R i-1 R i+1 • • • R M   .
begin the R i the i-th column of the matrix R. Notice now that:

det (R i ) = M j=2 (-1) j+1 ρ 1,j det(R i,j ), where det(R i,j ) is the (i, j)-th minor of R (since R is symmetric). Finally, det R = det (R) - M i,j=2 (-1) i+j ρ 1,i ρ 1,j det (R i,j ) = det (R)   1 - M i,j=2 (-1) i+j ρ 1,i ρ 1,j det (R i,j ) det (R)   = det (R)   1 - M i,j=2 (-1) i+j R 1,i R 1,j det (R i,j ) det (R) E 1 √ E i E j   = det (R) E 1   E 1 - M i,j=2 (-1) i+j R 1,i R -1 i,j R 1,j   = det (R) E 1 E 1 -r T R -1 r
So both criteria are equivalent if det(R) and E 1 are well defined and different than zero. What does it mean? Explain that correctly.

B Deriving formulaes for ∇ J and H J

In this section we will derive the formulae for the gradient and the Hessian of the criterion in [START_REF] Carbonetto | MATLAB primal-dual interior-point solver for convex programs with constraints[END_REF]. To do that we will need some matrix calculus rules. Let Y : R ⇒ R M ×M , be a matrix function depending on y, the following formulas hold:

• ∂ det (Y(y)) ∂y = det(Y(y))trace Y(y) -1 ∂Y(y) ∂y • ∂ trace (Y(y)) ∂y = trace ∂ Y(y) ∂y • ∂ Y(y) -1 ∂y = -Y(y) -1 ∂Y(y) ∂y Y(y) -1
Recall that the function we want to derivate is J = det R . We can then compute:

∂ J ∂t 1,k = ∂ det R ∂t 1,k = det R trace R-1 ∂ R ∂t 1,k . (24)
We can also compute the second derivative of the criterion:

∂ 2 J ∂t 1,j ∂t 1,k = ∂ ∂t 1,j det R trace R-1 ∂ R ∂t 1,k = det R trace R-1 ∂ R ∂t 1,j trace R-1 ∂ R ∂t 1,k + det R trace - R-1 ∂ R ∂t 1,j R-1 ∂ R ∂t 1,j + R-1 ∂ 2 R ∂t 1,j ∂t 1,k .
Hence, to be able to evaluate the gradient and the Hessian of J we need to compute the first and second derivatives of the matrix R. We will first compute these derivatives in function of the derivatives of R and r to finally compute the derivatives of R and r.

We can rewrite R as:

R = D     E 1 r T r R     D, (25)
where

D =     E -1/2 1 • • • 0 . . . . . . . . . 0 • • • E -1/2 M     . (26)
Notice that the matrix D does not depend on t. Hence the expressions for the derivatives of R are:

∂ R ∂t 1,k = D         0 ∂r ∂t 1,k T ∂r ∂t 1,k ∂R ∂t 1,k         D, (27)
INRIA and

∂ 2 R ∂t 1,j ∂t 1,k = D         0 ∂ 2 r ∂t 1,j ∂t 1,k T ∂ 2 r ∂t 1,j ∂t 1,k ∂ 2 R ∂t 1,j ∂t 1,k         D. (28)
So we just need the first and second derivatives of the vector r and the matrix R. Recall their expressions:

R =       R 2,2 (0) R 2,3 (t 1,3 -t 1,2 ) • • • R 2,M (t 1,M -t 1,2 ) R 2,3 (t 1,3 -t 1,2 ) R 3,3 (0) • • • R 3,M (t 1,M -t 1,3 ) . . . . . . . . . . . . R 2,M (t 1,M -t 1,2 ) R 3,M (t 1,M -t 1,3 ) • • • R M,M (0)       , and r =       R 1,2 (-t 1,2 ) R 1,3 (-t 1,3 ) . . . R 1,M (-t 1,M )      
.

The partial derivative of R is matrix filled with zeros expect for its k-th row and its k-th column that are equal to the following vector:

∂R ∂t 1k k =             R ′ 2,k (t 1,k -t 1,2 ) . . . R ′ k-1,k (t 1,k -t 1,k-1 ) 0 -R ′ k+1,k (t 1,k+1 -t 1,k ) . . . -R ′ M,k (t 1,M -t 1,k )             . ( 29
)
The partial derivative of r is:

∂r ∂t 1,k =         0 . . . -R ′ 1,k (-t 1,k ) . . . 0         (30)
We will differentiate two cases when computing the second derivative: j = k This will fill the diagonal of the hessian matrix. Notice that:

∂ 2 r ∂t 2 1,k =         0 . . . R ′′ 1,k (-t 1,k ) . . . 0         (31)
and that the partial second derivative of R is matrix filled with zeros expect for its k-th row and its k-th column that are equal to the following vector:

∂ 2 R ∂t 2 1k k =             R ′′ 2,k (t 1,k -t 1,2 ) . . . R ′′ k-1,k (t 1,k -t 1,k-1 ) 0 R ′′ k+1,k (t 1,k+1 -t 1,k ) . . . R ′′ M,k (t 1,M -t 1,k )             . ( 32
)
j > k This fills the lower triangular matrix of the hessian (and the upper triangular since we assume that the hessian is symmetric, i.e., J is twice continuously differentiable). The second derivative of r is null in this case, however the second derivative of R is not. Actually just two positions in the second derivative are not necessarily null: the jk-th and the kj-th being:

∂ 2 R ∂t 1,j ∂t 1,k kj = -R ′′ k,j (t 1,k -t 1,j ) (33)

C Estimating of the cross-correlation function

The previous derivation needs the cross-correlation function to be twice continuously differentiable. However, we just have the discrete signals. A natural question arises: how can we estimate the cross-correlation function from two discrete signals? We have chosen to compute the cross-correlation function of the interpolated signals. The following prsents how to estimate the cross-correlation function from some polynomial interpolation of two discrete signals.

Let x and y be two continuous signals. We'll denote their discretizations by x[n] and y[n]. We assume regular sampling at a rate F (and denote the sampling period T = 1/F ). The n-th sampling time will be denoted by t n = nT . We want to estimate the cross-correlation function of the signals x and y from the interpolation of x[n] and y[n].

INRIA

Let us assume a D-degree polynomial interpolation for each time interval:

x(t) = x n (t) = D d=0 a n,d (t -t n ) d ŷ(t) = y n (t) = D d=0 b n,d (t -t n ) d . ( 34
)
The task is to compute R xŷ (τ ). Notice that:

R xŷ (τ ) = x(t)ŷ(t -τ )dt (35)
= n tn+1 tn x(t)ŷ(t -τ )dt (36)
= n tn+1 tn x n (t)ŷ(t -τ )dt. (37)
Let us now decompose the delay in τ = mT + τ , with 0 ≤ τ < T . We have the following relations: We will know develop the two terms of the sum. The first term:

t ∈ (t n , t n + τ ) ⇒ t -τ ∈ (t n-m -τ , t n-m ) ⊂ (t n-m-1 , t n-m ) t ∈ (t n + τ , t n+1 ) ⇒ t -τ ∈ (t n-m , t n+1 -t m -τ ) ⊂ (t n-m , t n+1-m ) (38)
I 1 = tn+τ tn x n (t)y n-m-1 (t -τ )dt (40)
= D p,q=0 a n,p b n-m-1,q tn+τ tn (t -t n ) p (t -τ -t n-m-1 ) q dt (41) = D p,q=0 a n,p b n-m-1,q tn+τ tn (t -t n ) p (t -t m -τ -t n + m +T ) q dt (42) = D p,q=0 a n,p b n-m-1,q τ 0 t p (t + T -τ ) q dt (43) = D p,q=0 a n,p b n-m-1,q q k=0 q k τ 0 t p t k (T -τ ) q-k dt (44) = D p,q=0 a n,p b n-m-1,q q k=0 q k (T -τ ) q-k τ k+p+1 k + p + 1 (45) = D p,q=0 a n,p b n-m-1,q K 1 (q, p, τ , T ) (46)
And the second term:

I 2 = tn+1 tn+τ x n (t)y n-m (t -τ )dt (47) = D p,q=0 a n,p b n-m,q tn+1 tn+τ (t -t n ) p (t -τ -t n-m ) q dt (48) = D p,q=0 a n,p b n-m,q tn+1 tn+τ (t -t n ) p (t -τ -t n ) q dt (49) = D p,q=0 a n,p b n-m,q T τ t p (t -τ ) q dt (50) = D p,q=0 a n,p b n-m,q q k=0 q k T τ t p t k (-τ ) q-k dt (51) = D p,q=0 a n,p b n-m,q q k=0 q k (-τ ) q-k T k+p+1 -τ k+p+1 k + p + 1 (52) = D p,q=0 a n,p b n-m,q K 2 (q, p, τ , T ) (53)
We can plug these two equations into the original cross-correlation equation and we get:

R xŷ (τ ) = n D p,q=0 a n,p b n-m-1,q K 1 (q, p, τ , T ) + + D p,q=0 a n,p b n-m,q K 2 (q, p, τ , T ) = D p,q=0 R ap,bq [m + 1] K 1 (q, p, τ , T ) + R ap,bq [m] K 2 (q, p, τ , T )
From this equation we can also compute the derivative of the cross-correlation:

R ′ xŷ (τ ) = D p,q=0 R ap,bq [m + 1] ∂K 1 (q, p, τ , T ) ∂ τ + R ap,bq [m] ∂K 2 (q, p, τ , T ) ∂ τ (54)
as well as the second derivative:

R ′′ xŷ (τ ) = D p,q=0 R ap,bq [m + 1] ∂ 2 K 1 (q, p, τ , T ) ∂ τ 2 + R ap,bq [m] ∂ 2 K 2 (q, p, τ , T ) ∂ τ 2 . (55)
The formulas for the partial derivatives are:

∂K 1 ∂ τ = q k=0 q k (T -τ ) q-k-1 τ k+p T - τ q + p + 1 k + p + 1 (56) INRIA ∂ 2 K 1 ∂ τ 2 = q k=0 q k (T -τ ) q-k-2 τ k+p-1 τ 2 (q + p + 1)(q + p) k + p + 1 -2τ T (q + p) + T 2 (k + p) (57) ∂K 2 ∂ τ = q k=0 q k (-τ ) q-k-1 (k -q) T k+p+1 -τ k+p+1 k + p + 1 + τ k+p+1 . ( 58
) ∂ 2 K 2 ∂ τ 2 = q k=0 q k (-τ ) q-k-2 (1 + k -q)(k -q) k + p + 1 (T k+p+1 -τ k+p+1 ) + (k -p -2q)τ k+p+1 (59)
D Deriving the ∇∆ and H∆

∇∆ = 2 A, B -M 1 J T A (B -M 1 ) + J T B A - -( A 2 -1)J T B (B -M 1 ) -B -M 1 2 J T A A
where J A = -2νM ‡ and J B = -2ν 2 M ‡ diag(t), being ν the sound speed. We can also compute the Hessian of ∆:

E Analytic cross-correlation function

In this section we compute the closed-form cross-correlation function for two different signals:

The exponentially decaied sinusoidal:

x(t) = exp -t/λ sin(F t).

We finally have all the ingredients to compute the cross-correlation function: The full spectrum signal we consider has the following expression:

y(t) = N n=1 a i sin(F i t + φ i ) (65)
where the a i are drawn from a standard uniform distribution, the F i correspond to divisors of the sampling frequency F i = iF/L, where F is the sampling frequency and L is the acquisition length in samples (notice that in this case N < L/2 since higher frequencies will not be present due to the Nyquist theorem).

Notice that, we do not want to compute the expected cross-correlation of this signal, but the correlation of the acquired signals as if they were deterministic. This computation aims to have the ground truth to test the criterion we proposed, it is different from deriving general properties of this signal.

As in the previous case we can write:

  (a) Good (non-anomalous) TDEs (b) TDE standard deviation

Figure 1 :

 1 Figure 1: (Best seen in color) TDE method comparison based on random speech fragments. The plots show (a) the percentage of good TDEs and (c) the standard deviation as a function of SNR in an anechoic setup (T60 = 0).

  (a) the percentage of non-anomalous TDE is plot and in Fig. 2(b) the standard deviation of such time delays estimation is shown. Notice how the performance of both method improves with the SNR. Also, the proposed method INRIA (a) Non-anomalous TDE (b) TDE Standard Deviation

Figure 2 :

 2 Figure 2: Evaluation of the TDE performance of gTDE method. The x-axis corresponds to the SNR value (dB), the color to the methods (blue for gtde, green for bypairs), and the line style to the reverberation level (solid-circle for T60 = 0 and dashed-square for T60 = 0.1 s). (a) shows the percentage of non-anomalous TDE and (b) the standard deviation of this estimates.

Figure 3 :

 3 Figure 3: Evaluation of the localization performance of gTDE method. The x-axis corresponds to the SNR value (dB), the color to the methods (blue for gtde, green for bypairs), and the line style to the reverberation level (solid-circle for T60 = 0 and dashed-square for T60 = 0.1 s). (a) shows the percentage of localization inliers and (b) the mean squared error of localization error.

  That allows us to split the integral in two: R xŷ (t) = n tn+τ tn x n (t)y n-m-1 (t -τ )dt I1 + tn+1 tn+τx n (t)y n-m (t -τ )dt I2[START_REF] Yang | Cramer-Rao bound and optimum sensor array for source localization from time differences of arrival[END_REF]

H∆ = 2

 2 J T A (B -M 1 ) + J T B A J T A (B -M 1 ) + J T B A T + + A, B -M 1 J T A J B + D + J T B J A --2(J T B (B -M 1 ))(J T A A) T + ( A 2 -1)(E + J T B J B ) + + 2(J T A A)(J T B (B -M 1 )) T + B -M 1 2 J T A J A where D = -2ν 2 Diag (M † A) E = -2ν 2 Diag (M † (B -M 1 )).

R 2 2ME. 2

 22 x1,x2 (τ ) = λ 4 exp -t/λ (cos(F t) -F λ sin(F t)) 1 + (λF ) )/λ cos(F (t 2 + τ -t 1 )) exp -t/λ 2M 2mThe full spectrum sinusoidal

R

  y1,y2 (τ ) = M m y 1 (t)y 2 (t -τ )dt = N i,j=1 a i a j M m cos(F i (t -t 1 ) + φ i ) cos(F j (t -t 2 -τ ) + φ j )dtSo we need to solve the following general integral:I = cos(αt + β) cos(γ(t + δ) + ǫ)dt = 1 2 cos((α -γ) t + β -γδ -ǫ)dt α + γ) t + β + γδ + ǫ)dt I BThose two integrals can be expressed as:α -γ) t + β -γδ -ǫ) α -γ α = γ t cos(β -γδ -ǫ) α = γ α + γ) t + β + γδ + ǫ) α + γ α + γ = 0 t cos(β + γδ + ǫ) α + γ = 0(67)We can now compute the cross-correlation function:R y1,y2 (τ ) = N i,j=1

  (t -t 1 ) + φ i ) cos(F j (t -t 2 -τ ) + φ j )dt

				M
	a i a j cos(F = m N M i,j=1 j =i a i a j m M cos(F + i=1 a 2 i m a i a j cos(F = N M -t1
		i,j=1 j =i	m-t1
				M -t1
	+	i=1	a 2 i	m-t1
					t=M -t1
					+
					t=m-t1
	+	sin((F t=M -t1	+
					t=m-t1
					t=M -t1
	+		a 2		+
		i=1			t=m-t1
	+	sin(2(F 2F i	t=M -t1 t=m-t1

i i (t -t 1 ) + φ i ) cos(F j (t -t 2 -τ ) + φ j )dt + i (t -t 1 ) + φ i ) cos(F i (t -t 2 -τ ) + φ i )dt cos(F i t + φ i ) cos(F j (t + t 1 -t 2 -τ ) + φ j )dt + cos(F i t + φ i ) cos(F i (t + t 1 -t 2 -τ ) + φ i )dt = N i,j=1 j =i a i a j sin((F i -F j ) t + φ i -F m (t 1 -t 2 -τ ) -φ j ) F i -F j i + F j ) t + φ i + F j (t 1 -t 2 -τ ) + φ j ) F i + F j i t cos(F i (t 1 -t 2 -τ )) i t + φ i ) + F i (t 1 -t 2 -τ ))

RR n°7988

The full spectrum sinusoidal:

We want to compute R z1,z2 (τ ), where z is either x or y and z i means that the signal is delaied by t i . To be more precise, given the original signal z, we write the two following "acquisitions" of the signal:

and the goal is to derive a closed-form solution for R z1,z2 (τ ).

E.1 Exponentially decaied sinusoidal

We write, then:

We just need to compute the following general integral:

The integral above is computed by applyint integration by parts twice: cos(F t) exp -t/λ dt = λ 1 + (λF ) 2 exp -t/λ (cos(F t) -λF sin(F t)).

(64) INRIA