Hermine Biermé

Lionel Moisan

Frédéric Richard

A TURNING-BAND METHOD FOR THE SIMULATION OF ANISOTROPIC FRACTIONAL BROWNIAN FIELDS

Keywords: 1991 Mathematics Subject Classication. 60G60, 60G22, 60-04 Turning bands, simulation, fractional Brownian eld, Gaussian eld, anisotropy, non-stationarity, stationary increments. 1 HERMINE BIERMÉ 1,2 , LIONEL MOISAN 1 ,

. The TBM reduces the problem of simulating a eld in two dimensions by combining independent processes simulated on oriented bands. In the AFBF context, the simulation elds are constructed by solving an integral equation arising from the application of the TBM to non-stationary anisotropic elds. This garantees the convergence of simulations as their precision is increased. The construction is followed by a theoretical study of the convergence rate. Another key feature of this work is the simulation of band processes. Using self-similarity properties, processes are simulated exactly on bands with a circulant embedding method, so that simulation errors are exclusively due to the eld approximation. Moreover, we design a Dynamic Programming algorithm that selects band orientations achieving the optimal trade-o between computational cost and precision. Finally, we conduct a numerical study showing that the approximation error does not signicantly depend on the regularity of the elds to be simulated, nor on their degree of anisotropy. Experiments also suggest that simulations preserve eld statistical properties.

Introduction

In this paper, we address the issue of simulating realizations of a generic class of Gaussian elds, known as Anisotropic Fractional Brownian Fields (AFBF) and introduced in [START_REF] Bonami | Anisotropic analysis of some Gaussian models[END_REF]. Date: October 9, 2012.

The simulation of AFBF is an open issue whose complexity is mainly due to both the non-stationarity and the anisotropy of the elds. In [START_REF] Stein | Fast and exact simulation of fractional Brownian surfaces[END_REF], Stein described a specic method for the simulation of (isotropic) FBF. This method is based on a representation of the FBF by a locally stationary isotropic Gaussian eld, which is simulated using circulant embedding matrix techniques developped in [START_REF] Wood | Simulation of stationary Gaussian proceses on [0, 1] d[END_REF]. This simulation is exact and ecient on a regular grid. However, no locally stationary representation is available for general anisotropic elds so that Stein's method cannot be extended to this situation.

More generic methods based on covariance matrix factorizations [START_REF] Brouste | On fractional gaussian random elds simulations[END_REF][START_REF] Chan | An eective method for simulating Gaussian random elds[END_REF] can theoretically be applied to the AFBF simulation issue, but their computational cost is prohibitive and covariance functions are not known explicitely in the general case. Other methods based on the discretization of a continous spectral representation of the eld were used for the simulation of FBF in [START_REF] Saupe | The science of fractal images, chapter Algorithms for random fractals[END_REF] and AFBF in [START_REF] Ayache | Identication and series decomposition of anisotropic Gaussian elds[END_REF][START_REF] Biermé | Estimation of anisotropic gaussian elds through Radon transform[END_REF]. However, due to truncation or periodization of the spectral representation, the statistical properties of the simulated eld does not exactly match those of the theoretical eld.

In this paper, we focus on another generic simulation method, called the turning-band method (TBM) [START_REF] Journel | Geostatistics for conditional simulation of ore bodies[END_REF][START_REF] Matheron | The intrinsic random functions their application[END_REF]. The TBM essentially reduces the problem of simulating a eld in several dimensions to the problem of simulating several processes in one dimension. Indeed, consider the problem of generating a realization of a target eld X on a discrete set G of points of R 2 . Choose n lines (called turning bands) passing through a given origin and denote by θ i the angle indicating the direction of the ith band. The TBM is based on a combination of n appropriate processes (Y i ) 1≤i≤n independently simulated on each predened band:

(3)

∀ x ∈ G, X n (x) = n i=1 λ i Y i (x • u(θ i )),
where the λ i 's are positive weights and u(θ) = (cos(θ), sin(θ)) is the unit vector with direction θ. There are two major issues raised by the TBM. The rst one consists of determining appropriate weights λ i and band processes Y i which ensure the convergence of the turning-band eld X n to the target eld X as n tends to innity. The second one concerns the simulation of the processes Y i on the non-equispaced points {x•u(θ i ), x ∈ G}.

The convergence issue has been extensively studied in the case when the target eld is stationary [START_REF] Brooker | Two-dimensional simulation by turning bands[END_REF][START_REF] Christakos | Stochastic simulation of spatially correlated geo-processes[END_REF][START_REF] Dietrich | A simple and ecient space domain implementation of the turning bands method[END_REF][START_REF] Gneiting | Comment on "a simple and ecient space domain implementation of the turning bands method" by C. Dietrich[END_REF][START_REF] Gneiting | Closed form solutions of the two-dimensional turning bands equation[END_REF][START_REF] Mantoglou | Digital simulation of multivariate two-and three-dimensional stochastic processes with a spectral turning bands method[END_REF][START_REF] Matheron | The intrinsic random functions their application[END_REF]. In this case, the convergence can be obtained using stationary band processes. Let Cθ be the covariance of a band process in the direction θ. Take orientations (θ i ) 1≤i≤n uniformly distributed over (-π/2, π/2), and set λ i = π n . Then, the covariance of X n at point x is C n (x) = π n n i=1 Cθ i (x • u(θ i )) and, as a Riemann sum, it converges to T ( C)(x) = π 2 -π 2 Cθ (x • u(θ))dθ as n tends to innity. Hence, for the turning-band eld X n to converge to a target eld X of covariance C X , it suces that C satises the condition T ( C) ≡ C X . In the special case when the target eld is isotropic, Cθ ≡ C0 does not depend on θ, and the previous condition reduces to T (C 0 ) ≡ C X . This integral equation was solved explicitely for many dierent types of covariance C X (Gaussian, Whittle-Matérn, Cauchy, etc.), making possible to apply the TBM to a wide range of stationary isotropic elds [START_REF] Brooker | Two-dimensional simulation by turning bands[END_REF][START_REF] Christakos | Stochastic simulation of spatially correlated geo-processes[END_REF][START_REF] Dietrich | A simple and ecient space domain implementation of the turning bands method[END_REF][START_REF] Gneiting | Comment on "a simple and ecient space domain implementation of the turning bands method" by C. Dietrich[END_REF][START_REF] Gneiting | Closed form solutions of the two-dimensional turning bands equation[END_REF][START_REF] Lantuéjoul | Geostatistical Simulation: Models and Algorithms[END_REF]. In the anisotropic case, the HERMINE BIERMÉ 1,2 , LIONEL MOISAN 1 , AND FRÉDÉRIC RICHARD 3 equation was solved from spectral representations of covariances by expressing spectral densities of processes as a function of the one of the target eld [START_REF] Mantoglou | Digital simulation of multivariate two-and three-dimensional stochastic processes with a spectral turning bands method[END_REF].

The TBM can also be adapted to the simulation of non-stationnary elds with stationary increments. In such a situation, the convergence of the turning-band eld to the target one is rather expressed in terms of variograms. It can be obtained by taking band processes with stationary increments and variograms ṽθ (in the direction θ) which satisfy the condition (4)

v X (x) = π 2 -π 2 ṽθ (y • u(θ))dθ,
where v X is the variogram of the target eld (see Equation ( 6) later). In [START_REF] Yin | New methods for simulation of fractional Brownian motion[END_REF], the previous integral equation was solved in the particular case when v X is the variogram of an (isotropic) FBF of order H. It was shown that variograms v θ do not depend on θ and are proportional to the variogram of a one-dimensional FBM of the same order H. In [START_REF] Emery | A turning bands program for conditional co-simulation of cross-correlated Gaussian random elds[END_REF], this result was extended to other non-stationary isotropic Gaussian elds with spline generalized covariance. Similar ideas can be found in [START_REF] Dimitrakopoulos | Conditional simulation of intrinsic random functions of order k[END_REF][START_REF] Mantoglou | The turning bands method for simulation of random elds using line generation by a spectral method[END_REF][START_REF] Pardo-Igúzquiza | IRFK2D: a computer program for simulating intrinsic random functions of order k[END_REF] about the simulation of the so-called intrinsic random elds of order k generalizing elds with stationary increments [START_REF] Chilès | Geostatistics: modeling spatial uncertainty[END_REF][START_REF] Matheron | The intrinsic random functions their application[END_REF]. However, none of these works directly addressed the issue of simulating non-stationary anisotropic elds. In [START_REF] Setas | Modeling anisotropic and fractal two-dimensional elds: a tool for image simulation[END_REF][START_REF] Wu | Analyses and simulation of anisotropic fractal surfaces[END_REF], some attempts were made for the simulation of such elds with a TBM, but they only cover a few special cases. Let us also emphasize that the simulated elds are not Gaussian, so that Gaussian realizations may only be approximated by averaging several independent realizations and applying a Central Limit Theorem.

One of the main originalities of this paper is the construction of appropriate turningband elds for the simulation of AFBF, which are themselves centered Gaussian random elds with stationary increments. This is done by solving Equation (4) when v X is dened by Equations ( 1) and ( 2). This construction is completed by an analysis of the simulation error from both theoretical and numerical viewpoints. This study brings new insights into the TBM simulation error, which had been mainly investigated in the stationary isotropic case [START_REF] Chilès | Géostatique des phénomènes non stationnaires[END_REF][START_REF] Gneiting | Closed form solutions of the two-dimensional turning bands equation[END_REF][START_REF] Mantoglou | The turning bands method for simulation of random elds using line generation by a spectral method[END_REF].

In the construction of turning-band elds for AFBF, we show that band processes are one dimensional FBMs with varying Hurst indices. Hence, the application of the TBM to AFBF directly leads to the issue of simulating these processes on non-uniformly spaced points. In the context of TBM application, Yin [START_REF] Yin | New methods for simulation of fractional Brownian motion[END_REF] simulated FBMs using an adaptation of the spectral method [START_REF] Mantoglou | The turning bands method for simulation of random elds using line generation by a spectral method[END_REF][START_REF] Pardo-Igúzquiza | The Fourier integral method: an ecient spectral method for simulation of random elds[END_REF][START_REF] Shinozuka | Digital simulation of random processes and its applications[END_REF]. Based on a discretization of the spectral density of the process covariance, this method induces periodization eects and is inaccurate for mainly two reasons: the simulated process is not Gaussian and its covariance function only approximates the target one. In [START_REF] Emery | A turning bands program for conditional co-simulation of cross-correlated Gaussian random elds[END_REF][START_REF] Emery | TBSIM: A computer program for conditional simulation of threedimensional Gaussian random elds via the turning bands method[END_REF], periodization eects induced by Fourier methods was overcome using a continous spectral method. This method is fast and can be applied with arbitrary simulation points. However, it does not produce Gaussian realizations neither. In another context, Perrin et al. [START_REF] Perrin | Fast and exact synthesis for 1-D fractional Brownian motion a nd fractional gaussian noises[END_REF] developed a circulant embedding method (see [START_REF] Dietrich | Fast and exact simulation of stationary gaussian processes through circulant embedding of the covariance matrix[END_REF][START_REF] Wood | Simulation of stationary Gaussian proceses on [0, 1] d[END_REF] and Section 3.1) for the simulation of FBF. This method is fast and exact but requires equispaced simulation points. However, as shown in this paper, the issue of simulating FBM on band points can be reexpressed on equispaced points using self-similarity properties of FBM, as soon as bands orientations are conveniently chosen.

Hence, we can apply the circulant embedding method to obtain exact simulations of FBM on band points but the cost of these simulations depends on band orientations, and is higher than the one of the continuous spectral method in [START_REF] Emery | A turning bands program for conditional co-simulation of cross-correlated Gaussian random elds[END_REF][START_REF] Emery | TBSIM: A computer program for conditional simulation of threedimensional Gaussian random elds via the turning bands method[END_REF]. To reduce the global computational cost, we thus propose a Dynamic Programming [START_REF] Bellman | The theory of dynamic programming[END_REF] algorithm that selects band orientations in an optimal way. Gaussian eld X is characterized by its covariance function: (y, z) → Cov(X(y), X(z)).

A eld X has stationary increments if the law governing the eld

X(• + z) -X(z) is the same as X(•) -X(0) for all z ∈ R d .
When the eld X is centered and with stationary increments, we have ( 5)

∀y, z ∈ R 2 , Cov(X(y), X(z)) = v X (y) + v X (z) -v X (y -z),
where v X is the so-called variogram of X dened as ( 6)

∀y ∈ R 2 , v X (y) = 1 2 E((X(y) -X(0)) 2 ).
Hence, if the eld X is also Gaussian, its law is characterized by its variogram [START_REF] Biermé | Central limit theorems and quadratic variations in terms of spectral density[END_REF].

In this work, we deal with anisotropic fractional Brownian elds, which are centered Gaussian elds with stationary increments, characterized by a variogram of the form [START_REF] Averbuch | A framework for discrete integral transformations I -the pseudo-polar Fourier transform[END_REF] with a spectral density dened as in Equation [START_REF] Averbuch | A framework for discrete integral transformations II -the 2D discrete Radon transform[END_REF]. When c ≡ C > 0 and h ≡ H ∈ (0, 1)

are both constant, the variogram satises (see Remark 1.1.13 of [START_REF] Lacaux | Real Harmonizable multifractional Levy eld[END_REF] for instance)

(7) v(x) = 1 2 R 2 e ix•ζ -1 2 C ζ -2H-d dξ = C π 1 2 Γ(H + 1/2)γ(H) 2Γ(H + 1) x 2H ,
where for all H ∈ (0, 1),

γ(H) = π HΓ(2H) sin(Hπ) .

It follows that such elds are isotropic, which means that their law is invariant under rotation. They are also H-self-similar. When the function c is not constant but h ≡ H remains constant, the eld remains self-similar of order H but becomes anisotropic. When h is also allowed to vary, the eld is not self-similar anymore but, setting

H = essinf θ∈S d-1 ;c(θ)>0
h(θ), one can still nd a continuous modication of X with H as critical Hölder exponent [START_REF] Bonami | Anisotropic analysis of some Gaussian models[END_REF]. The fractal dimension of its graph is still linked with H by the relation D = 2 -H a.s. (see [START_REF] Xiao | Sample path properties of anisotropic Gaussian random elds[END_REF] for instance).

In general, it is dicult to get an explicit form of the AFBF variogram similar to the one expressed for the FBF in Equation [START_REF] Biermé | Estimation of anisotropic gaussian elds through Radon transform[END_REF]. However, we have computed explicitely the variogram of a particular class of AFBF which is slightly more general than the FBF. This eld, that we call elementary eld, is dened by a spectral density of the form [START_REF] Averbuch | A framework for discrete integral transformations II -the 2D discrete Radon transform[END_REF]

with c = 1 [α 1 ,α 2 ] for -π/2 ≤ α 1 < α 2 ≤ π/2,
and h ≡ H for H ∈ (0, 1). When α 2 = -α 1 = π/2, an elementary eld corresponds to a FBF of order H. As explained in Section 4, elementary elds will be of particular interest for the numerical evaluation of simulations.

Proposition 2.1. Let H ∈ (0, 1) and

-π/2 ≤ α 1 < α 2 ≤ π/2. Let denote v H,α 1 ,α 2 the variogram of an AFBF with h = H and c = 1 [α 1 ,α 2 ] . Then, ( 9
) ∀x ∈ R 2 , v H,α 1 ,α 2 (x) = 2 2H-1 γ(H)C H,α 1 ,α 2 (arg(x)) x 2H ,
where x is the Euclidean norm of x and C H,α 1 ,α 2 is a π-periodic function dened on

(-π/2, π/2] by C H,α 1 ,α 2 (θ) =            β H 1-sin(α 2 -θ) 2 + β H 1-sin(α 1 -θ) 2 if α 1 ≤ θ + π/2 ≤ α 2 β H 1+sin(α 2 -θ) 2 + β H 1+sin(α 1 -θ) 2 if α 1 ≤ θ -π/2 ≤ α 2 β H 1-sin(α 2 -θ) 2 -β H 1-sin(α 1 -θ) 2
otherwise with β H the Beta incomplete function given by

∀t ∈ [0, 1], β H (t) = t 0 u H-1/2 (1 -u) H-1/2 du,
and γ(H) is dened in Equation [START_REF] Biermé | ESAIM Proceedings: Mathematical Methods for Imaging and Inverse Problems, volume 26, chapter Anisotropic texture modeling and applications to medical image analysis[END_REF].

The proof of this proposition is given in appendix A. Now, let us consider the general case from which the TBM will follow.

v X (x) = 1 2 π/2 -π/2 γ(h(θ))c(θ)|x • u(θ)| 2h(θ) dθ,
where u(θ) = (cos(θ), sin(θ)) and γ(H) is dened in Equation [START_REF] Biermé | ESAIM Proceedings: Mathematical Methods for Imaging and Inverse Problems, volume 26, chapter Anisotropic texture modeling and applications to medical image analysis[END_REF].

Proof. Let x ∈ R 2 . Then acording to (1) and ( 2),

2v X (x) = R 2 e ix•ζ -1 2 c(arg(ζ))|ζ| -2h(arg(ζ))-2 dζ = 2π 0 +∞ 0 e ir(x•u(θ)) -1 2 c(θ)r -2h(θ)-1 drdθ,
by a change of variables in polar coordinates. But, for H ∈ (0, 1) and t ∈ R,

+∞ 0 e irt -1 2 r -2H-1 dr = 1 2 R e ist -1 2 |s| -2H-1 ds = 1 2 γ(H)|t| 2H ,
according to (7.2.13) of [START_REF] Samorodnitsky | Stable non-Gaussian random processes[END_REF]. Then the result follows by π-periodicity of h and c.

The integral equation ( 10) is of the form (4) with ṽθ ( θ) . This means that ṽθ is a solution of the integral equation (4) when v X is the variogram of an AFBF. Now recall that a FBM of order H is a centered Gaussian process with stationary increments and variogram w H (t) = 1 2 |t| 2H for all t ∈ R. Hence, ignoring the factor γ(h(θ))c(θ) depending on the orientation θ, the variogram ṽθ is equal to the one of a FBM of order h(θ), also varying with θ.

•) = γ(h(θ))c(θ) 1 2 | • | 2h(
According to previous remarks, we now specify turning-band elds for AFBF simulations. Given an ordered set Θ = (θ i ) 1≤i≤n of band orientations -π/2 ≤ θ 1 < . . . < θ n ≤ π/2, and a set Λ = (λ i ) 1≤i≤n ∈ [0, +∞) n of appropriate band weights, turning-band elds have the form [START_REF] Brouste | On fractional gaussian random elds simulations[END_REF]

X Θ,Λ (x) = n i=1 λ i γ(h(θ i ))c(θ i )Y i (x • u(θ i )), ∀x ∈ R 2 ,
where the Y i 's are n independent FBM of order h(θ i ). In the remaining of the text, turning-band elds X Θ,Λ will be called the simulation elds and processes Y i will always be FBF of order h(θ i ). We will also describe the precision of simulation elds X Θ,Λ using the variable ( 12)

ε Θ = max i=1,••• ,n+1 (θ i -θ i-1 ) , with n = |Θ|,
and θ 0 ∈ [-π/2, θ 1 ] and θ n+1 ∈ [θ n , π/2] are xed directions chosen according to the AFBF function c. Remark that a uniform choice for the orientations consists in choosing

θ n = θ n+1 = π/2 and θ i = -π/2 + iπ/n for 0 ≤ i ≤ n. This choice leads to ε Θ = π n while, for general orientations, we always have ε Θ ≥ θn-θ 0 n ≥ θn-θ 1 n .
The error of simulating X by X Θ,Λ may be expressed, at point x ∈ R 2 , as the Kolmogorov distance between the random variables X Θ,Λ (x) and X(x), that is, [START_REF] Chilès | Geostatistics: modeling spatial uncertainty[END_REF]

d Kol (X Θ,Λ (x), X(x)) = sup t∈R |P(X Θ,Λ (x) ≤ t) -P(X(x) ≤ t)| .
When this distance tends to 0, it implies that the random variable X Θ,Λ (x) tends to X(x) in distribution. As stated next, due to our Gaussian framework, this distance can be further bounded by a distance between variograms of simulation and target elds at x. Theorem 2.3. Let X Θ,Λ be a simulation eld dened as in Equation [START_REF] Brouste | On fractional gaussian random elds simulations[END_REF]. Then, X Θ,Λ is a centered Gaussian random eld on R 2 with stationary increments and variogram [START_REF] Chilès | Géostatique des phénomènes non stationnaires[END_REF]

v Θ,Λ (x) := n i=1 λ i γ(h(θ i ))c(θ i )w h(θ i ) (x • u(θ i )).
Moreover, X Θ,Λ (0) = X(0) = 0 a.s. and, for all x = 0, [START_REF] Christakos | Stochastic simulation of spatially correlated geo-processes[END_REF]

d Kol (X Θ,Λ (x), X(x)) ≤ 2 |v X (x) -v Θ,Λ (x)| v X (x) . Choosing (Θ n , Λ n ) n in such a way that v Θn,Λn (x) -→ n→+∞ v X (x) for all x ∈ R 2 , (X Θn,Λn (x)) x∈R 2 f dd -→ n→+∞ (X(x)) x∈R 2 ,
where f dd -→ stands for convergence of nite dimensional distributions.

Proof. Let us write

X θ i (x) := Y i (x•u(θ i )), for x ∈ R 2 , with Y i a FBM of order h(θ i ). First, remark that X θ i (0) = Y i (0) = 0 a.s. Moreover, since Y i is a centered Gaussian random HERMINE BIERMÉ 1,2 , LIONEL MOISAN 1 , AND FRÉDÉRIC RICHARD 3 process it is clear that X θ i is a centered Gaussian random eld on R 2 . Finally, since Y i has stationary increments, for any x 0 ∈ R 2 , writing t 0,i = x 0 • u(θ i ) ∈ R, {X θ i (x + x 0 ) -X θ i (x 0 ); x ∈ R 2 } = {Y i (x • u(θ i ) + t 0,i ) -Y i (t 0,i ); x ∈ R 2 } f dd = {Y i (x • u(θ i )) -Y i (0); x ∈ R 2 }, = {X θ i (x) -X θ i (0); x ∈ R 2 }.
It follows that X Θ,Λ is a centered Gaussian random eld on R 2 with stationary increments as a sum of independent centered Gaussian random elds on R 2 with stationary increments elds. Since X Θ,Λ (0) = 0 a.s.,

v X Θ,Λ (x) = 1 2 E X Θ,Λ (x) 2 = 1 2 Var(X Θ,Λ (x)), since X Θ,Λ is centered, = 1 2 n i=1 λ i γ(h(θ i ))c(θ i )Var(Y i (x • u(θ i ))), by independence, = 1 2 n i=1 λ i γ(h(θ i ))c(θ i )|x • u(θ i )| 2h(θ i ) = v Θ,Λ (x).
Let N ∼ N (0, 1), then for x = 0,

d Kol (X Θ,Λ (x), X(x)) = d Kol ( v Θ,Λ (x)N, v X (x)N ).
Then [START_REF] Christakos | Stochastic simulation of spatially correlated geo-processes[END_REF] follows from the fact that d Kol (σ N, σN ) ≤ 2 |σ-σ | σ , already remarked in [START_REF] Biermé | Central limit theorems and quadratic variations in terms of spectral density[END_REF]. Actually, there is nothing to prove when |σ-σ | σ > 1 2 . Otherwise we use the fact that for σ > 1 and z > 0 one has

P(z < N ≤ σz) ≤ (σ -1)ze -z 2 /2 ≤ σ 2 -1. Now, let assume that (Θ n , Λ n ) is such that v Θn,Λn (x) -→ n→+∞ v X (x) for all x ∈ R 2 .

By stationarity of the increments, for all

n ≥ 1, for all x, y ∈ R 2 , Cov(X Θn,Λn (x), X Θn,Λn (y)) = v Θn,Λn (x) + v Θn,Λn (y) -v Θn,Λn (x -y),
and similarly for Cov(X(x), X(y)) and v X . It follows that Cov(X Θn,Λn (x), X Θn,Λn (y)) tends to Cov(X(x), X(y)) for all x, y ∈ R 2 . Using a Cramér-Wold device, this implies that the eld (X Θn,Λn (x)) x∈R 2 converges to (X(x)) x∈R 2 , for nite dimensional distributions.

Let us quote that, since v Θ,Λ appears as a numerical approximation of the integral giving v X , one can choose (Θ n , Λ n ) n in such a way that v Θn,Λn (x) tends to v(x) for x ∈ R 2 . This implies that d Kol (X Θn,Λn (x), X(x)) → 0 so that X Θn,Λn (x) tends to X(x) in distribution.

Note that conversely, since X Θn,Λn (x) and X(x) are centered Gaussian variables, v Θn,Λn (x) tends to v(x) as soon as X Θn,Λn (x) tends to X(x) in distribution. The next section is devoted to the rate of convergence. 2.3. Approximation error. We can choose Θ, Λ such that the following uniform bounds hold for approximation of elementary elds. Proposition 2.4. Let c and h be two π-periodic mesurable functions dened on

(-π/2, π/2] by h = H for some H ∈ (0, 1) and c = 1 [α 1 ,α 2 ] for -π/2 ≤ α 1 < α 2 ≤ π/2. Let Θ = (θ i ) 1≤i≤n with α 1 ≤ θ 1 < . . . < θ n ≤ α 2 and θ 0 = α 1 , θ n+1 = α 2 . Choose Λ as (16) λ 1 = θ 2 -θ 0 and λ i = θ i+1 -θ i for 2 ≤ i ≤ n.
Then, one can nd a positive constant C > 0, independent of Θ, Λ, such that for all

x ∈ R d , (17)
d Kol (X Θ,Λ (x), X H,α 1 ,α 2 (x)) ≤ Cε min(2H,1) Θ ,
where the precision parameter ε Θ is dened in Equation [START_REF] Chan | An eective method for simulating Gaussian random elds[END_REF].

Moreover, when choosing ( 18)

λ 1 = (θ 1 -θ 0 )+ θ 2 -θ 1 2 , λ n = (θ n+1 -θ n )+ θ n -θ n-1 2
and λ i = θ i+1 -θ i-1 2 for 2 ≤ i ≤ n-1, one can nd a positive constant C > 0, independent of Θ, Λ, such that for all x ∈ R d , ( 19
) d Kol (X Θ,Λ (x), X H,α 1 ,α 2 (x)) ≤ C    ε 3 Θ δ -2+min(2H,1) Θ + ε 1+min(2H,1) Θ if H = 1/2, ε 3 Θ δ -1 Θ | log(δ Θ )| + ε 2 Θ if H = 1/2 with δ Θ = min 1≤i≤n-1 (θ i+1 -θ i ).
The proof of Proposition 2.4 is postponed to Appendix B. In [START_REF] Emery | A spectral approach to simulating intrinsic random elds with power and spline generalized covariances[END_REF], the authors also propose a TBM to synthesize isotropic FBF (case

α 1 = -π/2 and α 2 = π/2) in general dimension d ≥ 2.
However, the processes simulated on the bands are not Gaussian so that the Kolmogorov distance between the simulated random variable and

X H,-π/2,π/2 (x)
is bounded by the Berry Esseen bound given by n -1/2 , with n the number of bands, (see Equation ( 27) in [START_REF] Emery | A spectral approach to simulating intrinsic random elds with power and spline generalized covariances[END_REF]). Moreover, in their case this distance also depends on the

point x ∈ R 2 .
Let us compare with our results. Note that when orientations are chosen

uniformly one has ε Θ = δ Θ = α 2 -α 1 n
so that, choosing a rectangular rule, we obtain in [START_REF] Dietrich | A simple and ecient space domain implementation of the turning bands method[END_REF] a bound given by n -min(2H,1) , while for a trapezoidal rule, we obtain in [START_REF] Dimitrakopoulos | Conditional simulation of intrinsic random functions of order k[END_REF] a bound given by n -1-min(2H,1) when H = 1/2 and n -2 log(n) when H = 1/2. Moreover let us emphasize that our bounds do not depend on x ∈ R 2 . This could be generalized to other self-similar AFBF (h ≡ H) under regularity assumptions on c. In the general case, our bounds depend on x through the term 1/v X (x). However, we obtain uniform bounds for the dierence |v X (x) -v Θ,Λ (x)| when x is in a compact set, as stated in the next proposition.

Proposition 2.5. Let assume that h and c are piecewise

C 1 on (-π/2, π/2]. Let Θ = (θ i ) 1≤i≤n with -π/2 ≤ θ 1 < . . . < θ n ≤ π/2
containing the singular points of h and c and θ 0 = -π/2, θ n+1 = π/2. Let T be a compact set of R d . Then, one can nd

Λ = (λ i ) 1≤i≤n ∈ [0, +∞) n and a positive constant C T > 0, independent of Θ, Λ, such that for all x ∈ T , ( 20
)
|v X (x) -v Θ,Λ (x)| ≤ C T ε min(2H,1) Θ , where H = min θ∈[-π/2,π/2]
h(θ) > 0 and ε Θ is dened in Equation [START_REF] Chan | An eective method for simulating Gaussian random elds[END_REF]. If moreover, h and

c are piecewise C 2 on (-π/2, π/2], one can nd Λ = (λ i ) 1≤i≤n ∈ [0, +∞) n and a positive constant C T > 0, independent of Θ, Λ, such that for all x ∈ T , ( 21
)
|v X (x) -v Θ,Λ (x)| ≤ C T    ε 3 Θ δ -2+min(2H,1) Θ + ε 1+min(2H,1) Θ if H = 1/2, ε 3 Θ δ -1 Θ | log(δ Θ )| + ε 2 Θ if H = 1/2 with δ Θ = min 1≤i≤n-1 (θ i+1 -θ i ).
The proof of Proposition 2.5 is postponed to Appendix B. A bound for the Kolmogorov distance is then obtained using the fact that

d Kol (X Θ,Λ (x), X(x)) ≤ 2 |v X (x)-v Θ,Λ (x)| v X (x)
. Then we may choose an increasing sequence of (Θ n ) n and choose (Λ n ) n as given in Proposition 2.5. Therefore, if ε Θn → 0, the sequence of random elds (X Θn,Λn (x)) x∈R 2 converges to

(X(x)) x∈R 2 ,
for nite dimensional distributions in view of Theorem 2.3. The next section is devoted to the simulation of the band processes.

3. Fast and exact simulation on bands 3.1. Simulation of 1D fractional Brownian motions. Several methods for the synthesis of 1D fractional Brownian motions have been proposed in the literature. Most of them are approximate procedures. However, considering equispaced points on the band one can get exact simulations using the circulant embedding method [START_REF] Dietrich | Fast and exact simulation of stationary gaussian processes through circulant embedding of the covariance matrix[END_REF]. Let us briey recall this procedure. Let H ∈ (0, 1) and B H a fractional Brownian motion. We consider

Z H = (B H (t + 1) -B H (t)
) t∈R , the fractional Gaussian noise that is a stationary process with covariance function given by

Cov(Z H (t), Z H (s)) = r H (|t -s|) with r H (t) = 1 2 (|t + 1| 2H -2|t| 2H + |t -1| 2H ).
Let l ≥ 1, then the vector (Z H (0), . . . , Z H (l -1)) is a centered Gaussian vector of size l with Toeplitz covariance matrix R H (l) = (r H (|i -j|)) 0≤i,j≤l-1 . One can embed R H (l) in a circulant matrix of size 2l given by S H (l) = circ(s H (l)) with s H (l) = (r H (0), . . . r H (l -

1) r H (l -2), . . . r H (1)).
The main interesting property of circulant matrices is that they are diagonalized in the discrete Fourier basis, with their eigenvalues given by the Discrete Fourier Transform of their rst row. In particular one has

S H (l) = 1 2l F * 2l diag(F 2l s H (l))F 2l , where F 2l = (e iπjk l ) 0≤i,j≤2l-1 .
The main results of [START_REF] Perrin | Fast and exact synthesis for 1-D fractional Brownian motion a nd fractional gaussian noises[END_REF] is that for all H ∈ (0, 1) and l ≥ 1, F 2l s H (l) always has positive entries, so that S H (l) is a covariance matrix. Moreover, if ε (1) 2l and ε (2) 2l are independent vectors of law N (0, I 2p ) the vectors

Z (1) = 1 √ 2l F * 2l diag(F 2l s H (l)) 1/2 (ε (1) 2l + iε (1) 2l
and

Z (2) = 1 √ 2l F * 2l diag(F 2l s H (l)) 1/2 (ε (1) 2l + iε (2) 2l
are independent with common law N (0, S H (l)). In particular, one has

(Z H (0), . . . , Z H (l -1)) d = (Z (i) 0 , . . . , Z (i)
l-1 ) for i = 1, 2 and using stationarity of the increments of B H and the fact that B H (0) = 0 a.s., one has, for all m ≤ l,

(B H (k)) -m≤k≤l-m = j<k+m Z H (j) - j<m Z H (j) -m≤k≤l-m
, with the convention that j<0 = 0. Let us emphasize that this procedure is very fast since, choosing l as a power of 2, the cost is reduced to O(l log(l)) using the Fast Fourier Transform algorithm.

3.2. Choice of bands. We consider the exact simulation of X Θ,Λ on the discrete grid r -1 Z 2 ∩ [0, 1] 2 for some r ≥ 1. Then for any i with 1 ≤ i ≤ n we have to simulate on each band of direction u(θ i )

Y i (x • u(θ i )); x ∈ r -1 Z 2 ∩ [0, 1] 2 = Y i k 1 r cos(θ i ) + k 2 r sin(θ i ) ; 0 ≤ k 1 , k 2 ≤ r .
Note that when θ i = π/2 we can simply use the fact that

Y i k 2 r ; 0 ≤ k 2 ≤ r d = r -h(θ i ) {Y i (k 2 ) ; 0 ≤ k 2 ≤
r}, by self-similarity. When cos(θ i ) = 0 we may choose θ i such that tan θ i = p i q i , with p i ∈ Z and q i ∈ N and use that fact that

Y i k 1 r cos(θ i ) + k 2 r sin(θ i ) ; 0 ≤ k 1 , k 2 ≤ r d = cos(θ i ) rq i h(θ i ) {Y i (k 1 q i + k 2 p i ) ; 0 ≤ k 1 , k 2 ≤ r} .
Thus, the band with direction u(θ i ) involves the simulation of a 1D fractionnal Brownian motion on a discrete interval of length r(|p i | + q i ). The computational cost of this simula-

tion is O(C(|p i | + q i ))
, where O(C(l)) is the computational cost of the Fourier Transform used in the procedure described in Section 3.1 to simulate a 1D fractional Brownian motion on the discrete interval {0, . . . , rl}. If the Fast Fourier Transform with powers of two is used, then one has C(l) = 2 log 2 (rl) log 2 (rl) , where x denotes the upper integer part of x. Finally, the overall simulation process has to nd a trade-o between the computational cost

C(Θ) = i C(|p i | + q i )
and the precision of the simulation, which is controlled by

E(Θ) = ε min(2H,1)
Θ . The optimal choice of Θ is discussed in the following section.

3.3. Band selection by Dynamic Programming. As we just saw in the previous section, we need to restrain our choice of band orientations to angles θ that correspond to vectors (q, p) with integer coordinates, that is, such that tan θ = p q . Moreover, in order to control the total computational cost, we would like to favor small factors (small values of |p| + q) while controlling the repartition of bands in order to keep ε Θ small. A simple (but non-optimal) solution to select the set of angles Θ = (θ i ) i consists in using a uniform discretization θi = α 1 + i n (α 2 -α 1 ), then choosing for each i a rational approximation p i q i of tan θi (this can be done very eciently using the appropriate convergent of the continued fraction associated to tan θi ). However, as we shall see now, one can nd a simple algorithm, based on Dynamic Programming [START_REF] Bellman | The theory of dynamic programming[END_REF], that is able to select a set of angles Θ that minimizes the computational cost C(Θ) under the error control constraint

ε Θ ≤ ε.
In practice, we restrain our choice to angles that can be written under the form θ = ∠(q, p), where (p, q) belongs to

V N = {(p, q); -N ≤ p ≤ N, 1 ≤ q ≤ N, gcd(p, q) = 1, α 1 < ∠(q, p) < α 2 } and ∠(q, p) is the measure in [-π 2 , π 2 ]
of the angle of the vector (q, p), obtained by ∠(q, p) = arctan p q . The integer N should be chosen large enough to ensure that the optimal solution only involves vectors from V N . It seems that choosing N = 1 + 1 tan ε (so that ∠(N, 1) < ε) is enough, though we do not have a proof of this (even if this were not true, the algorithm we present here would yield slightly sub-optimal sets of bands, with little practical consequences). Now assume that the set V N has been sorted into a sequence (p k , q k ) 1≤k≤n such that the associated angular sequence

θ k = ∠(q k , p k ) is increasing. Writing e k = C(r(|p k | + q k ))
the elementary cost associated to a band with orientation θ k , we can rewrite the total computational cost of a set of angles Θ = (θ i k ) 1≤k≤s as

C(Θ) = s k=1 e i k .
We add the convention that θ n+1 = α 2 and θ 0 = α 1 (with the associated elementary cost e 0 = 0). Now, for 0 ≤ i ≤ n + 1, let us call c i the minimal cost that can be realized with a sequence i 1 = i, i 2 , . . . i s = n + 1 for some integer s. Then, c 0 is the optimal cost we look for, and for all 0 ≤ i ≤ n we have [START_REF] Emery | A spectral approach to simulating intrinsic random elds with power and spline generalized covariances[END_REF] c i = e i + min j; j>i, θ j ≤θ i +ε c j .

This induction formula (called Bellman Equation in the framework of Dynamic Programming) permits us to compute the optimal costs c n , c n-1 , . . . c 0 recursively (the initialization being made with c n+1 = 0). Moreover, each time the minimum in ( 22) is computed, we consider one optimal index

k i ∈ arg min j; j>i, θ j ≤θ i +ε c j ,
then an optimal sequence i 1 , i 2 , . . . i s can be computed by tracking back indexes that achieve the optimal cost c 0 . This sequence is given by

i 1 = k 0 , i 2 = k i 1 , . . . i s = k i s-1 ,
where the value of s is obtained using the fact that i s+1 = n + 1. In the end, the desired sequence of integer vectors is simply (p k , qk ) 1≤k≤s , where (p k , qk ) = (p i k , q i k ) for all

1 ≤ k ≤ s.
The whole procedure we just described is given in pseudo-code in Algorithm 1 of Appendix C.

Numerical Study

This section is devoted to the numerical evaluation of anisotropic fractional Brownian eld (AFBF) simulations obtained by turning bands.

Let us recall some notations. The eld X is the theoretical eld to be simulated (AFBF). Its variogram v X is of the form (1) with a spectral density dened by [START_REF] Averbuch | A framework for discrete integral transformations II -the 2D discrete Radon transform[END_REF]. The eld X Θ,Λ is the turning-band simulation eld dened by Equation ( 11) for some sets Θ and Λ giving band orientations and weights, respectively. The variogram v Θ,Λ of X Θ,Λ is dened by Equation ( 14).

In all experiments, the set Θ of band orientations was computed automatically using the Dynamic Programming algorithm described in Section C with a constraint on eld precision. The precision parameter ε Θ associated to the set Θ is dened as in Equation [START_REF] Chan | An eective method for simulating Gaussian random elds[END_REF]. Weights λ i of Λ are dened to fulll the condition (16) of Proposition 2.4. 4.1. The use of elementary elds. Our evaluation was focused on elementary elds whose spectral density is given by Equation (2) taking h ≡ H for some H ∈ (0, 1) and

c = 1 [-α,α] for some 0 < α ≤ π/2.
Elementary elds are specied by only two parameters, H and α, which can be interpreted as regularity and anisotropy parameters, respectively. The Hölder regularity of these elds being equal to H (see [START_REF] Bonami | Anisotropic analysis of some Gaussian models[END_REF] for instance), it increases as H tends to 1. When α = π/2, elementary elds correspond to usual isotropic fractional Brownian elds of Hurst index H. When 0 < α < π/2, these elds are no longer isotropic. In this case, they are some kind of fractional Brownian elds whose non-null frequency components are restricted between frequency directions -α and α. As α decreases to 0, non-null eld frequency components become more and more focused around the horizontal direction.

On Figure 2, some elementary eld realizations are shown for illustrating both the eect of increasing H on the eld regularity and the eect of decreasing α on its anisotropy. For the evaluation, we considered elementary elds of varying degrees of regularity and anisotropy, taking all parameter pairs (H, α) for H in {0.2, 0.5, 0.8} and α in {π/6, π/3, π/2}.

Let us further mention that variograms of elementary elds can be computed using the closed form given in Equation ( 9). As it will appear next, this is of particular interest for the computation of evaluation criteria. On Figure 3, some of these variograms are presented for dierent degrees of regularity and anisotropy. Finally, let us notice that any anisotropic fractional Brownian eld whose spectral density is dened with piecewise constant functions h and c can be decomposed as a sum of independent elementary elds. Hence, although achieved on elementary elds, our evaluation accounts for more general anisotropic fractional Brownian elds.

4.2. Approximation error. As mentionned in Section 2.3, simulation errors result from the distance separating simulation and theoretical elds. This distance can be dened as the Kolmogorov distance between distributions of theoretical and simulation elds at each position x. As stated in Theorem 2.3 (Equation ( 15)), the Kolmogorov distance

is further bounded by d Θ,Λ (x) = |v X (x)-v Θ,Λ (x)| v X (x)
, which is proportional to the error made when approximating the variogram of X by the one of X Θ,Λ . Moreover, when X is an elementary eld, it is possible to compute the bound d Θ,Λ (x) using closed forms of v Θ,Λ (x)

and v X (x) given by Equations ( 14) and ( 9), respectively. Hence, using elementary elds, we could numerically evaluate a simulation error by averaging values of d Θ,Λ (x) over points

x of a uniform subgrid of [0, 1] 2 :

(23)

d Θ,Λ = p k,l=1 d Θ,Λ k p , l p = p k,l=1 |v X ( k p , l p ) -v Θ,Λ ( k p , l p )| v X ( k p , l p )
, with p = 64.

As this is evidenced by Equation ( 17) of Proposition 2.4, the measured error d Θ,Λ depends on the precision parameter ε Θ of the simulation eld. Figure 4 illustrates the eect of increasing ε Θ (i.e reducing the simulation precision) on simulations of a fractional Brownian eld of Hurst index H = 0.2. When the precision becomes too low (ε Θ ≥ 0.25), eld realizations have some stripes in dierent directions, and simulation eld variograms present some singularities on lines radiating from the origin. This well-known eect, often called artifact banding in the literature [START_REF] Gneiting | Closed form solutions of the two-dimensional turning bands equation[END_REF][START_REF] Mantoglou | The turning bands method for simulation of random elds using line generation by a spectral method[END_REF][START_REF] Emery | TBSIM: A computer program for conditional simulation of threedimensional Gaussian random elds via the turning bands method[END_REF][START_REF] Emery | A turning bands program for conditional co-simulation of cross-correlated Gaussian random elds[END_REF], is due to the fact that the contribution of a band process Y i to the sum dening the simulation eld (see Equation [START_REF] Brouste | On fractional gaussian random elds simulations[END_REF]) is null for points on the line orthogonal to the band direction θ i and passing through origin.

On Figure 5, we plotted values of error bounds d Θ,Λ obtained for dierent elementary elds as a function of the precision parameter ε Θ . Whatever the eld parameters, error bounds varied almost linearly with respect to ε Θ . They did not seem to depend on the regularity parameter H. However, they were slightly dependent on the anisotropy parameter α, especially at low precision (ε Θ > 0.03). Error bounds of all elds fell below 1% when ε Θ < 0.02, such a precision being reached with around 150 simulation bands.

4.3. Estimation error. We also conducted numerical experiments to evaluate errors which arise when estimating eld features (e.g. parameters, variograms) from eld simulations.

Applying the TBM (with 1321 bands on a 64 × 64 grid of [0, 1] 2 ), we simulated 2000

independent realizations {y (k) , k = 1, • • • , 2000} of a given elementary eld X of variogram v X . Given k ∈ {1, • • • , 2000}
, we then computed the empirical variogram v (K) (x) at position x using the K rst samples: (1)

v (K) (x) = 1 2K K k=1 (y (k) (x)) 2 .
Figure 7. Variations of the estimation error relative to the number K of samples for elementary elds with dierent values of H and α: (1)

H = 0.2, (2)
H = 0.5, (3)
H = 0.8.
random variables, and H x 0 be the hypothesis that their distribution is the same as the one of X(x), i.e. a centered gaussian distribution with variance 2v X (x). For testing H x 0 , we dene the rejection interval {D (K) (x) > c}, where the statistic

D (K) (x) = |V (K) (x)-v X (x)| v X (x) with V (K) (x) = 1 2K K k=1 (Y (k) (x)) 2 . Under assumption H x 0 , K V (K) (x) v X (x) has a χ 2 distribution
of degree K. Hence, under H x 0 , the probability of the rejection interval can be computed, and the rejection bound c can be set according to a level of test. Given a realization

d (k) (x) of D (K) (x)
, it also possible to evaluate the p-value p (K) (x) of the test (i.e the minimal risk of rejecting H x 0 ) as

p (K) (x) = 1 -P (Z (K) < K(d (k) (x) + 1)) + P (Z (K) < K(-d (k) (x) + 1)),
where Z (K) is a χ 2 K random variable. To evaluate the simulation-to-model adequacy, we computed for dierent K the average p (K) of p-values p (K) (x) at positions x of the 64 × 64 grid of [0, 1] 2 . On Figure 8, the mean p-values are plotted as a function of the sample number K. Whatever the value of parameters H and α, mean p-values are all above 0.3, indicating that hypotheses H x 0 of adequacy are not rejected at low risks. Besides, the mean p-values seems to reach an upper bound which is below 1 (around 0.5). This is probably due to both the approximation error and estimator inacurracy. (1)

Discussion

We have constructed turning-band elds suited to the simulation of AFBF. This construction is based on the resolution of an integral equation specic to the non-stationary anisotropic context of AFBF. This ensures the convergence of simulation elds to target elds as the precision increases. Moreover, the band processes involved in the denition of simulation elds are simulated exactly using a circulant embedding method. Hence, errors produced by the simulation method are exclusively due to the approximation of the target eld by the simulation eld. This approximation error was evaluated theoretically and numerically. From a numerical point of view, we observed that it does not depend signicantly on the regularity of target elds, nor on their degree of anisotropy. Experiments have also suggested that simulations preserve the statistical properties of the target eld. Besides, we obtained good simulation results with few bands (around 150) at a low computational cost.

The evaluation was achieved on some elementary elds. However, simulation possibilities oered by the TBM go far beyond those elds, as there is a large choice of parameter denitions and tunings. Using the TBM, it chiey becomes possible to visualize truthfully realizations of dierent anisotropic eld models studied in the literature [START_REF] Bonami | Anisotropic analysis of some Gaussian models[END_REF][START_REF] Biermé | Central limit theorems and quadratic variations in terms of spectral density[END_REF][START_REF] Davies | Fractal analysis of surface roughness by using spatial data[END_REF][START_REF] Richard | Statistical tests of anisotropy for fractional Brownian textures. application to full-eld digital mammography[END_REF]. In the generic model dened by Equation ( 2), we recall that the eld anisotropy is introduced through two direction-dependent and π-periodic functions: the Hurst index function h and the topothesy function c. So as to illustrate the eect of varying these parameter functions, we considered three functions of increasing regularity:

• a discontinuous function g 1 : µ 1 , µ 2 ∈ (0, 1), g 1 (ω) = µ 1 if ω ∈ (-π 4 , π 4 ) and g 1 (ω) = µ 2 if ω ∈ (-π, π)\(-π 4 , π 4 ),
• a continuous but not dierentiable function g 2 : for µ 1 , µ 2 ∈ (0, 1), g 2 (0) = µ 1 , g 2 (-π/2) = g(π/2) = µ 2 , and g 2 is piecewise linear on (-π/2, 0) and (0, π/2).

• an innitely dierentiable function Fixing the Hurst index function to a constant (h ≡ 0.5), we also simulated eld realizations with topothesy functions c = g i for i = 1, 2, 3 and dierent pairs of parameter values (µ 1 , µ 2 ); results are shown on Figure 10. In these realizations, the degree of anisotropy can be dened as the dierence µ 2 -µ 1 between maximal and minimal topothesy values.

From a column to the next one, it is increased, while the eld regularity remains the same for all realizations (H = 0.2). As it can be observed by comparing realizations on a same row, variations of the anisotropy degree cannot be visually detected on textures.

However, comparing realizations on a same column, we can notice that the regularity of the topothesy function has an eect on the eld texture. As previously, some line patterns are present on textures when the topothesy function is discontinous.

In the simulations we presented, eld realizations were generated on a regular subgrid of [0, 1] 2 . Using our TBM approach, it is however possible to simulate elds on other sets of non-uniformly spaced positions. To do so, the only condition is that position coordinates and (µ 1 , µ 2 ) = (0.2, 0.8), respectively. are all rational; this is required for the exact simulation of fractional Brownian motions on turning bands (refer to Section 3.1). The pseudo-polar grid is an example of a set of points satisfying this simulation condition [START_REF] Averbuch | A framework for discrete integral transformations I -the pseudo-polar Fourier transform[END_REF]. Such a grid is of particular interest for computing discrete Radon transforms [START_REF] Averbuch | A framework for discrete integral transformations II -the 2D discrete Radon transform[END_REF], as its points are uniformly spread on dierent lines radiating from the origin. But Radon transforms are one of the key features for the construction of parameter estimators for AFBF [START_REF] Biermé | Estimation of anisotropic gaussian elds through Radon transform[END_REF][START_REF] Richard | Statistical tests of anisotropy for fractional Brownian textures. application to full-eld digital mammography[END_REF]. Hence, those estimators could be better discretized and evaluated using simulations on a pseudo-polar grid. Due to the ability of the TBM to simulate elds on quasi-arbitrary points, it also becomes possible to simulate eld deformations. For instance, let A be the 2 × 2-matrix of an ane transform (with rational components) and X = X • A the deformation of the random eld X by the ane transform A. Realization of X on a uniform grid G of [0, 1] 2 can be obtained by applying the TBM to the simulation of X on {Ax, x ∈ G}. Figure 11 presents an illustration corresponding to the deformation of a fractional Brownian eld of Simulation is a central issue concerning the investigation of anisotropic elds. As shown previously, a simulation technique such as the TBM can serve as a tool for visualizing mathematical properties of anisotropic models under study. From an application viewpoint, it can also help assessing the similarity between model realizations and real-world images. Besides, having reliable simulations is critical for the evaluation of model parameter estimators. In future works, we plan to use TBM simulations of AFBF to evaluate the estimators we constructed using quadratic variations [START_REF] Biermé | Estimation of anisotropic gaussian elds through Radon transform[END_REF][START_REF] Richard | Statistical tests of anisotropy for fractional Brownian textures. application to full-eld digital mammography[END_REF]. We also intend to use those simulations to rene the adequacy between models and radiographic images we analyze for the characterization of osteoporosis and breast cancer [START_REF] Biermé | Parametric estimation for gaussian operator scaling random e lds and anisotropy analysis of bone radiograph textures[END_REF][START_REF] Biermé | ESAIM Proceedings: Mathematical Methods for Imaging and Inverse Problems, volume 26, chapter Anisotropic texture modeling and applications to medical image analysis[END_REF][START_REF] Richard | Statistical tests of anisotropy for fractional Brownian textures. application to full-eld digital mammography[END_REF].

θ i+1 θ i (g x (θ) -g x (θ i )) dθ ≤ c 2 (θ i+1 -θ i ) 1+min (1,2H) .

Then, choosing λ 1 = θ 2 -θ 0 and λ i = θ i+1 -θ i for 2 ≤ i ≤ n, one has

θ i+1 θ i g x (θ) - g x (θ i ) + g x (θ i+1 ) 2 dθ ≤ sup θ∈[θ i ,θ i+1 ] |g x (θ)| (θ i+1 -θ i ) 3 12 . (30)
Note also that using [START_REF] Lacaux | Real Harmonizable multifractional Levy eld[END_REF] one always has [START_REF] Mantoglou | The turning bands method for simulation of random elds using line generation by a spectral method[END_REF] θ i+1 θ i g x (θ) -g x (θ i ) + g x (θ i+1 ) 2 dθ ≤ c 2 (θ i+1 -θ i ) 1+min (1,2H) . which gives the result using ( 30) with ( 32) and ( 33) for the two rst sums and using [START_REF] Lantuéjoul | Geostatistical Simulation: Models and Algorithms[END_REF] and [START_REF] Mantoglou | The turning bands method for simulation of random elds using line generation by a spectral method[END_REF] for the other terms. The general case where arg(x) = π/2 can be computed similarly.

It follows that

Proof of Proposition 2.5. The proof is similar to the proof of Proposition 2. instead of g x (θ) given by [START_REF] Kolmogorov | Wienersche Spiralen und einige andere interessante Kurven in Hilbertsche Raum[END_REF]. Note that when h and c are assumed of class C l (l = 1 or h(θ). These estimates allow to proceed as in the proof of Proposition 2.4. The result follows by summing the integrals over which the functions h and c are regular using the fact that H ≤ H(α 1 , α 2 ) for all α 1 , α 2 .

Figure 1 .

 1 Figure 1. Comparison of an approximate uniform sampling of turning bands and their selection by dynamic programming.

Figure 2 .

 2 Figure 2. Realizations on [0, 1] 2 of elementary elds obtained for dierent values of H and α using the TBM with 5900 bands.

Figure 3 .

 3 Figure 3. Variograms on [0, 1] 2 of elementary elds for dierent values of H and α.

Figure 4 .

 4 Figure 4. Eect of reducing the eld precision on simulations of a fractional Brownian eld of Hurst index H = 0.2. On the rst row, realizations of simulation elds X Θ,Λ of decreasing precision and, on the second row, corresponding variograms: (1) ε Θ = 0.04 (n = 103), (2) ε Θ = 0.077 (n = 51), (3) ε Θ = 0.25 (n = 15) and (4) ε Θ = 0.464 (n = 7).

Figure 6 .

 6 Figure 6. Eect of increasing the number of samples on simulations of a fractional Brownian eld of Hurst index H = 0.2 (simulations were done with 1321 bands on a 64 × 64 grid of [0, 1] 2 ): Empirical variograms v (K) computed with (1) K = 60, (2) K = 1000, (3) K = 2000 samples, and (4) theoretical variogram.

Figure 8 .

 8 Figure 8. Mean p-values p (K) of adequacy tests relative to the number K of samples for elementary elds with dierent values of H and α: (1) H = 0.2, (2) H = 0.5, (3) H = 0.8.

g 3 :

 3 g 3 (ω) = µ 1 r(ω) + µ 2 (1 -r(ω)) with r(ω) = (1 + sin(2ω + π/2))/2 for ω ∈ [0, π/2], and g 3 (ω) = g 3 (-ω) for ω ∈ [-π/2, 0].

Fixing the topothesy function

  to a constant (c ≡ 1), we rst simulated eld realizations with Hurst index functions h = g i for i = 1, 2, 3 and dierent pairs of parameter values (µ 1 , µ 2 ); results are shown on Figure9. In these realizations, the degree of anisotropy can be measured as the dierence µ 2 -µ 1 between maximal and minimal Hurst indices. It is the same for realizations of the rst and second columns (µ 2 -µ 1 = 0.3), and higher for those of the third column (µ 2 -µ 1 = 0.5). Moreover, the Hölder regularity of those realizations is equal to µ 1 . It is the same for realizations of the rst and third columns (µ = 0.2) and higher for those of the second column (µ = 0.5). Comparing realizations on a same row, we clearly see the eect of anisotropy and regularity variations on eld textures: as the eld regularity decreases, the texture gets rougher, and, as the eld anisotropy increases, texture patterns get more obviously oriented. Besides, comparing realizations on a same column, we can observe texture dierences induced by changing the regularity of the Hurst index function h in the model. In particular, realizations obtained with a discontinous function h (on the rst row) have some linear patterns which are not present on those obtained with a more regular function h (on the second and third rows).

Figure 9 .

 9 Figure 9. Field realizations obtained with dierent Hurst index functionsh. For all realizations, the topothesy function c ≡ 1. On the rst, second and third rows, Hurst index functions are h = g 1 (discontinuous), h = g 2 (continuous but not dierentiable), and h = g 3 (innitely dierentiable), respectively. On the rst, second and third columns, Hurst index functions are specied by parameter pair values (µ 1 , µ 2 ) = (0.2, 0.5), (µ 1 , µ 2 ) = (0.5, 0.8),

Figure 10 .

 10 Figure 10. Field realizations obtained with dierent topothesy functions. For all realizations, the Hurst index function h ≡ 0.2. On the rst, second and third rows, topothesy functions are c = g 1 (discontinuous), c = g 2 (continuous but not dierentiable), and c = g 3 (innitely dierentiable), respectively. On the rst, second and third columns, topothesy functions are specied by parameter pair values (µ 1 , µ 2 ) = (1, 5), (µ 1 , µ 2 ) = (1, 100), and (µ 1 , µ 2 ) = (1, 1000), respectively.

  example also shows that deformating isotropic elds is a means to construct anisotropic elds.

Figure 11 .

 11 Figure 11. Shear of a fractional Brownian eld of Hurst index 0.5.

3 First 2 α 1 | 1 . 1 α 2 α 1 |

 3211121 let us remark that since θ → | cos(θ-θ )| 2H is continuous, non negative and non identically equal to 0 one can nd c 1 > 0 such that for allx ∈ R 2 , α cos(θ -arg(x))| 2H dθ ≥ c It follows that d Kol (X Θ,Λ (x), X H,α 1 ,α 2 (x)) ≤ c -1 cos(θ -arg(x))| 2H dθ -n i=1 λ i | cos(θ i -arg(x))| 2H .Now, let us write (26)g x (θ) = | cos(θ -arg(x))| 2H ,and remark that (27)|g x (θ) -g x (θ )| ≤ 2|θ -θ | min(2H,1) for all θ, θ ∈ R,using the fact that ||t| 2H -|t | 2H | ≤ 2|t -t | min(2H,1) for all t, t ∈ [-1, 1].It follows that there exists c 2 > 0 such that for 0 ≤ i ≤ n,

α 2 α 1 g

 1 x (θ) -n i=1 λ i g x (θ i ) ≤ c 2 (α 2 -α 1 ) max 0≤i≤n (θ i+1 -θ i ) min(1,2H) . Moreover g x is of class C 2 on R {arg(x) -π/2 + πZ} with (29) |g x (θ)| ≤ c 3 | cos(θ -arg(x))| 2H-2 , for all θ / ∈ arg(x) -π/2 + πZ,for some c 3 > 0 (non depending on x). According to the trapezoidal rule, when [θ i , θ i+1 ] ∩ {arg(x) -π/2 + πZ} = ∅,

  -θ i ) g x (θ i ) + g x (θ i+1 ) 2 + (θ 1 -θ 0 )g x (θ 1 ) + (θ n+1 -θ n )g x (θ n ) = x (θ) -g x (θ i ) + g x (θ i+1 ) 2 dθ + θ 1 θ 0 g x (θ) -g x (θ 1 )dθ + θ n+1 θn g x (θ) -g x (θ n )dθ ≤ x (θ) -g x (θ 1 )dθ + θ n+1 θn g x (θ) -g x (θ n )dθ ,

  4, consideringgx (θ) = γ(h(θ))c(θ)|x • u(θ))| 2h(θ) = γ(h(θ))c(θ) x 2h(θ) | cos(θ -arg(x))| 2h(θ) .

2) on [α 1 , α 2 ]

 12 ⊂ [-π/2, π/2], one can nd C > 0 such that for all x ∈ T , | gx (θ) -gx (θ )| ≤ C|θ -θ | min(2H(α 1 ,α 2 ),1) , for all θ, θ ∈ [α 1 , α 2 ],and when l = 2,| gx (θ)| ≤ C| cos(θ -arg(x))| 2H(α 1 ,α 2 )-2 , for all θ ∈ [α 1 , α 2 ] with θ / ∈ arg(x) -π/2 + πZ, where H(α 1 , α 2 ) = min θ∈[α 1 ,α 2 ]

c n+1 ← 0 6 for 7 cmin ← ∞ 8 j 10 if c j ≤ cmin then 11 cmin ← c j 12 k 20 while i ≤ n do 21 i ← k i 22 (

 678101112202122 i = n, n -1, . . . 0 do ← i + 1 9 while j ≤ n + 1 and θ j ≤ θ i + ε do p s , qs ) ← (p i , q i ) 23 s ← s + 1 24 end 25 return (p k , qk ) 1≤k≤s26

MAP5, UMR-CNRS 8145, Université Paris Descartes, PRES Sorbonne Paris Cité, 45 rue des Saints-Pères, 75006 Paris, France.E-mail address: [bierme,moisan]@mi.parisdescartes.fr

LMPT, UMR-CNRS 7350, Fédération Denis Poisson FR-CNRS 2964, Université François- Rabelais de Tours, Parc de Grandmont

, 37200 Tours, France.3 LATP, UMR-CNRS 7353, Aix-Marseille Université, 39 rue F. Joliot Curie, 13453 Marseille Cedex 13, France. E-mail address: frederic.richard@univ-amu.fr

Acknowledgements

This work is part of the research program MATAIM, supported by the Agence Nationale pour la Recherche (ANR-09-BLAN-0029-01).

Finally, we computed an estimation error

at point x, and its average over points x of the grid (24)

, with p = 64.

On Figure 6, some empirical variograms v (K) estimated from simulations of a fractional Brownian eld of Hurst index H = 0.2 are compared to the theoretical variogram v X of X. This illustrates both the convergence of empirical variograms to the theoretical one as K tends to +∞, and estimation errors due to the lack of samples.

On Figure 7, we plotted the estimation error d (K) as a function of the sample number K for dierent elementary elds. For a same value of the regularity parameter H, the convergence of the error to zero is about the same for all values of the anisotropy parameter α. However, the convergence gets slower and slower as H increases. In all cases, around 1000 samples are required for the error to get below 5%.

Besides, we built a statistical test to check the adequacy of simulations to the model.

For any position

We will use the following lemma.

Proof. Since we assume that

which gives the result.

This allows us to get the next result, which concludes the proof.

according to Lemma A.1. Note that by a change of variables, for all θ ∈ R ( 25)

and the result follows.

Let assume thatπ 2 + 2Zπ = ∅. Then, one can nd k ∈ Z such that a + 2kπ < -π/2 < b + 2kπ, and according to Lemma

which concludes this case using [START_REF] Journel | Geostatistics for conditional simulation of ore bodies[END_REF]. The last case is similar.

Appendix B. Proofs of section 2.3

Proof of Proposition 2.4. Note that X Θ,Λ (0) = X H,α 1 ,α 2 (0) = 0 a.s. so that d Kol (X Θ,Λ (0), X H,α 1 ,α 2 (0)) = 0. Let x ∈ R 2 with x = 0. Then, the error of approximation is bounded by Build the set V of all integer vectors (q, p) ∈ {1, . . . , N } × {-N, . . . , N } such 2 that gcd(p, q) = 1 and α 1 < ∠(q, p) < α 2 Sort V into a sequence (p k , q k ) 1≤k≤n with k → θ k := ∠(q k , k k ) increasing 3

Compute the angle θ k := ∠(q k , p k ) and the cost e k associated to each (q k , p k )

Add extremal angles: θ 0 ← α 1 (e 0 ← 0) and θ n+1 ← α 2


Urtzi Etxeberria
email: u.etxeberria@gmail.com

Ricardo Etxepare
email: retxepare@gmail.com

When quantifiers do not agree: three systems

Keywords: Weak quantifiers, agreement alternation, plurality

Basque weak quantifiers optionally agree with the inflected verb in number. This papers' main aim is to study the dialectal variation shown by this phenomenon. The study will show that it is necessary to differentiate at least three systems: the western-central system, one that we will call the transition system, and the eastern system (souletin). The western-central system allows the presence of non-agreeing weak quantifiers in every case-marked position, ergative, dative or absolutive; the transition system does not allow it with ergative case arguments, and the oriental system allows it only with absolutive case arguments. In the latter system, the distribution of non-agreeing quantifiers is identical to that of bare nouns: bare nouns are only possible in those positions where absolutive case is assigned.

Introduction: the phenomenon

Basque plurality denoting nominal expressions trigger obligatory agreement in number with the inflected verb. i

(1) Anek liburu-ak erosi ditu/*du Ane.erg book-D-pl buy have.pl/have.sg 'Ane has bought (the) books' This rule has an exception in so called 'vague' weak quantifiers in Basque, which optionally agree in number with the inflected verb (2a-d) (see [START_REF] Rotaetxe | Numeral eta aditzen arteko konkordantziaz[END_REF][START_REF] Txillardegi | Oinarri bila[END_REF][START_REF] Txillardegi | Oinarri bila[END_REF]EGLU 1985;[START_REF] Etxepare | Euskal Gramatika Lehen Urratsak I[END_REF]. ii (2) a. Bezero asko etortzen da/dira halako egunetan customer many come-hab is/are such days-in 'A lot of customers come in such days' b. Bezero gehiegik eskatu du/dute arrain zopa customer too-many-erg asked aux-sg/aux-pl fish soup 'Too many customers asked for fish soup' c. Maiak lagun gutxi ikusi du/ditu gaur Maia.erg friend few seen aux.sg/aux.pl today 'Maia has seen few students today'

The notion of what we mean by 'vague' weak quantifier can be intuitively grasped by means of the following contrast:

(3) a. Mila ikasle etorri dira/*da thousand student come aux-pl/aux-sg 'One thousand students came' b. Milaka ikasle etorri dira/da thousand-suffix student come aux-pl/sg 'Thousands and thousands of students came' Whereas (3a), which involves a 'definite' quantity, triggers plural agreement in the inflected verb, (3b), which involves a non-definite quantity (equivalent to thousands of in English), only optionally triggers agreement. Cardinal quantifiers always trigger plural agreement in standard basque. Vague quantificational expressions constructed out of them, on the other hand, may not.

This phenomenon is general in the Basque area, with some interesting and systematic dialectal variation that we will try to synthesize here. The present paper offers a descriptive account of the variation involved in optional number agreement in the Basque area, as well as some basic generalizations that provide syntactic cues for a unified analysis. A full syntactic explanation of the dialectal variation related to this phenomenon is beyond the scope of this paper. One solid conclusion that follows from our discussion is that non-agreeing quantificational expressions are not counting expressions, but rather expressions related to what Borer (2005) has called a "stuff divider": a functional head whose semantic contribution is to portion out the denotation of count terms so that they can interact compositionally with the counting function. In that context, vague quantifiers merely measure the noun. Measures constitute the other quantificational domain in Basque that presents an agreement alternation in number. iii (4) Hiru litro ardo edan du/ditu three liter wine drunk aux-sg/aux-pl 'He/she drank three liters of wine'

We may wonder at this point what the agreement alternation is: is it an alternation between plural number features and singular ones? Or is the singular agreement form just a default, selected in the absence of any number feature? It is not easy to answer to this question by looking at the inflected forms directly. However, if we move to other syntactic contexts, the answer seems to favor the conclusion that third singular agreement, in the context of vague quantifiers in Basque, is just a default, with no correspondence with actual number features.

One such context is provided by secondary predication, which requires agreement in number (see Artiagoitia, 1994). The example in (5) gives an illustrative example with a Small Clause complement.

( Whereas a vague quantifier that agrees in plural with the inflected verb licenses a secondary predicate with a plural suffix -k on it, a vague quantifier that does not agree in plural can not license singular agreement in the secondary predicate either. The conclusion seems to be that agreement in singular with the quantifiers that do not agree in plural with the verb is impossible, and that therefore, the relevant quantifier forms must lack number features, either plural or singular. iv That the problem is in number agreement and not, say, in the ability of non-agreeing quantifiers to license a secondary predication is shown by the following fact: if we allow for a secondary predicate that does not have number, secondary predication with vague quantifiers becomes possible. One relevant configuration involves the [-ta] suffix, an adverbial ending that attaches to participles, which does not agree in number in Basque. When the participial substitutes for the [determiner+number] suffix, secondary predication with vague quantifiers becomes possible (7).

(7) Liburu asko hondatu-ta ikusi dut/ditut book many worn-out.part seen aux.sg/aux.pl 'I've seen many books worn-out'

The paper is organized as follows: In section 2 we present the received analysis concerning the agreement alternation in Basque. Section 3 provides arguments against this view. Sections 4-7 show the properties of Basque non-agreeing quantifiers and their dialectal variation. We distinguish three systems: (i) central, (ii) transition system (Lapurdian), (iii) eastern (Souletin). In section 8 we show the similarities between the non-use of the article in both Souletin and in some Romance languages. This section allows us to state a general syntactic condition on the non-agreeing cases. Section 9 concludes the paper.

A previous view: non-agreeing quantifiers are masses

The descriptive grammar of Euskaltzaindia (1985: 223-224) assimilates the absence of number agreement with weak quantifiers to the absence of number in mass terms. Take for instance the contrast in (8). According to the Academy's grammar, (9a) and (9b) do not have the same interpretation:

whereas "in the first case we consider a mass of books; in the other case we consider one book and then another one, and another one, and so on" (1985: 223). To make things clearer, the grammar presents the following case. In (10a) harri 'stone' is taken to be non-count, as a big quantity of stone. In (10b) it refers to a big quantity of stones (as a count term). The Academy's grammar does not go beyond the intuition above. Although we will not pursue this line of analysis, we share the intuition that (10b) offers more opportunities for an individualized treatment of the stone than (10a). For instance, (10b) would be more appropriate to describe a situation where demonstrators attack the police by throwing stones to them. This implies the existence of individualized pieces of stone, and a multiplicity of stone-throwing events. (10a) on the other hand, would be more appropriate to describe loads of stone being dumped during some road construction. For argument's sake, if we were to reformulate the Academy's proposal slightly, it could be stated as saying that number morphology coerces masses into counts (11), whereas absence of number morphology coerces count nouns into masses (12).

Are non agreeing quantifiers mass?

It can be shown however that non-agreeing quantifiers are not mass terms. As a starting point, we consider Pelletier's well known thought experiment (1975) to characterize mass terms. He proposes the existence of two imaginary machines, that he calls the Universal Grinder and the Universal Objectifier. For the Universal Grinder, we are to imagine a device which can grind anything, no matter how big or small. Into one end of the device "is inserted an object of which some count expression is true, and from the other end spews forth the finely-ground matter of which it is composed. So a hat is entered into the grinder and after a few minutes there is hat all over the floor" (from Pelletier and Schubert 1989:342). This is so despite the fact that we could also have said that there is felt all over the floor, using a mass expression. Examples of this type "show that many count expressions can be seen to already have within them a mass sense or a mass use" (ibidem: 343). Taking the word sagar 'apple' as our putative count term, we could take (13) to involve the mass coming out of the Universal Grinder.

(13) Entsaladak sagar pixkat dauka salad-D-erg apple bit has 'The salad has a bit of apple in it' Take, however, something like ( 14), with a non-agreeing vague quantifier. The sentence in ( 14), with a non-agreeing quantifier, does not involve a mass term, in Pelletier's sense: what I have seen in ( 14) is not scattered pieces of student, but a number of students, all of them of a piece. True, the force of this argument against a mass-approach to non-agreeing quantifiers depends on the force of Pelletier's metaphor to characterize mass terms as a whole. We know that in this sense, the metaphor is not comprehensive enough.

Other mass terms appear to reflect objects that we would better locate in the entering side of the machine. This is the case of mass terms like furniture or crockery (Chierchia 1998):

ground-up furniture and furniture do not mean the same, despite the mass status of the term.

In any case, even with simple ambiguous nouns such as apple, the mass-approach falls short of accounting for the range of interpretations that non-agreeing cases have. Consider a sentence like ( 15).

(15) Plater honetan sagar asko ikusten dut dish this-in apple many see aux-sg (i) 'I see a lot of apple in this dish'

(ii) 'I see a lot of apples in this dish'

As shown by the translations, non-agreeing quantifiers can be interpreted in two ways:

either as mass terms, referring to a quantity of apple, or as referring to a plural set of (whole) apples. In other words: the sentence in (15) can be interpreted as making reference to, say, a dish containing a set of piled-up entire apples. The mass-approach has nothing to say about this second interpretation.

Other properties distinguishing mass terms from non-agreeing cases lead us to reject the mass approach to non-agreeing quantifiers. Lonning (1987) Finally, we note that some of the quantifiers that give rise to the alternation just cannot quantify over mass terms. This is the case of zenbait 'some' and hainbat 'a sizeable quantity'.

(17) shows that even the non-agreeing cases do not support a mass interpretation. √ 'We brought some meats'

Up until now, we have concentrated on showing the differences that exist between nonagreeing quantifiers and mass terms. In the sections that follow, we will mainly concentrate on the dialectal variation that non-agreeing quantifiers show, and making as thorough a description as possible of this variation. As will be made clear, there are at least three systems in Basque when it comes to the distribution of non-agreeing quantifiers: Central-western, Transitional (Lapurdian), and Eastern (Souletin).

The distributive nature of non-agreeing quantifiers

One of the characterizing properties of non-agreeing quantifiers (which further distinguishes them from mass terms) is their distributive nature [START_REF] Etxepare | Euskal Gramatika Lehen Urratsak I[END_REF]. They can only be interpreted distributively, and this sets certain restrictions on the kind of predicate they can attach to.

Distributive readings

Consider for instance the contrast between ( 19) and ( 20).

(

Transition system: Lapurdian

Syntactic distribution of non-agreeing quantifiers

The transition system shows some differences compared to the Central-western system when it comes to the distribution of non-agreeing quantifiers. In the Central-western system non-agreeing quantifiers are grammatical in all grammatical functions, whereas in the Transition system this is not so: non-agreeing quantifiers can appear in S position, but only with absolutive case (29a), they don't accept to appear with the ergative case (29b); they can appear in IO position, with dative case (29c); and they can also appear in DO position with absolutive case (29d). Thus, non-agreeing quantifiers appear to be unable to appear with the ergative case.

Reciprocals

Reciprocals come to show exactly the same thing, that is, non-agreeing quantifiers have a distributive nature in this system. As was the case in the Central-western system, nonagreeing quantifiers are incompatible with reciprocals, as the ungrammaticality of (35b)

shows.

Enumeration and anaphora

As expected, in the transition system non-agreeing quantifiers cannot make reference to specific individuals. As a consequence, non-agreeing quantifiers cannot be antecedent to anaphors (36b) and they don't allow the enumeration of individuals, i.e. it is possible to make reference to the members of the set we are talking about, (37b).

(36) a. Bezero anitz i sartu dira gaur. _ i ez dira oso pozik atera customer many come aux.pl today neg aux.pl very happy leave Thus, the main difference between the Central-western and the Transition system is related to the possibility of non-agreeing quantifiers to appear with the ergative case. Non-agreeing quantifiers in the former system have no problem to appear with ergative case whereas in the latter system, non-agreeing quantifiers cannot take ergative case. 6. Eastern system: "Souletin"

Syntactic distribution of non-agreeing quantifiers

For the third system, we follow the description provided by [START_REF] Coyos | Le parler basque souletin des Arbailles: Une approche de l'ergativité[END_REF] for the dialect of Arbailles. In this system, as was the case in the Transition system, non-agreeing quantifiers cannot appear in all grammatical functions. They can appear in S position, but only with absolutive case (39a), not with ergative case (38a). They cannot appear in IO position, with dative case (40). A quote by Coyos is in order here: "N avec le datif et déterminé par un quantificateur indéfini: si l'indice de datif est present dans le syntagme verbal, ce sera celui avec le pluriel". Finally, non-agreeing quantifiers have no problem to appear in DO position;

in fact, the non-agreeing case is more commonly used than the agreeing one (41

Rest of dialects

The rest of Basque dialects do not accept BNs and the presence of the article is necessary if the sentence is going to be grammatical (cf. Laka 1993, Artiagoitia 1997[START_REF] Artiagoitia | The functional structure of Basque noun phrases[END_REF]. viii The presence of this overt number marker is closely related to the presence of the definite article [-a] (cf. [START_REF] Etxeberria | Quantification and domain restriction in Basque[END_REF]Etxeberria , 2011, in prep), in prep). [START_REF] Etxeberria | Izen eta Gertakarien Gaineko Kuantifikazioa[END_REF], 2009, in

Conclusions

Basque weak quantifiers display a number agreement alternation with the inflected verb.

This paper has investigated the dialectal variation of this phenomenon in Basque. We have seen that at least three different systems must be distinguished: central-western, transitional (Lapurdian), and eastern (Souletin). The paper has centered on the following grammatical issues and their geographical distribution: (i) syntactic contexts where the absence of agreement is allowed: in the central-western variety non-agreeing quantifiers can appear in all grammatical functions; the transitional system does not allow non-agreeing quantifiers with ergative case; and the eastern system only allows non-agreeing quantifiers in DO position; and (ii) the parallel distribution of non-agreeing quantifiers and bare nouns in Souletin:

Souletin, and only Souletin, allows BNs in Basque, and the syntactic distribution of these BNs is parallel to the one shown by non-agreeing cardinal quantifiers.

  ) Liburuak hondatu(*-ak) ikusi ditut book-D.pl worn-out.pl seen aux.pl 'I've seen (the) books worn-out' The sentence (5) contains a Small Clause predicate hondatuak 'worn-out' which obligatorily agrees in number with the subject liburuak 'books'. Now consider the contrast in (6). (6) a. Liburu asko hondatuak ikusi ditut book many worn-out.pl seen aux.pl 'I've seen many books worn-out' b. *Liburu asko hondatua ikusi dut book many worn-out.sg seen aux.sg 'I've seen many books worn-out'

( 8

 8 ) a. Haragi asko jaten du meat much eat-hab aux.sg 'He eats a lot of meat' b. Haragi asko jaten ditu meat many eat-hab aux.pl 'He eats many types of meat'The presence of number agreement in (8b) triggers a count interpretation of the mass term haragi 'meat', which comes to denote a set of individualized meat types. The grammar of the Academy suggests that the absence of number agreement with count terms has the opposite effect: it converts count terms into mass terms. The grammar comments on the following sentences in (9).

( 14 )

 14 Ikasle asko ikusi dut gaurko batzarrean student a lot of seen I-have today's meeting-D-in 'I have seen a lot of students in today's meeting'

  (17) a. Zenbait ardo edan dugu some wine drunk aux-sg * 'We drank some wine' √ 'We drank some wines' b. Hainbat haragi ekarri dugu some meat brought aux-sg * 'We brought some meat'

  4. Non-agreeing quantifiers: Central-western system4.1. Syntactic distribution of non-agreeing quantifiersIn this system non-agreeing quantifiers can occur in all syntactic positions and in all grammatical functions: in subject (S) position, both with ergative or absolutive case (18a); in indirect object (IO) position, with dative case (18b); and in direct object (DO) position with absolutive case. (18) a. Subject: Erg.: Azkenean gazte asko-k altxatu behar izan zuen harria finally young a lot of-erg lift must have aux-sg stone-D 'Ultimately, many youngsters had to lift the stone' Abs.: Ikasle asko etorri da festara student a lot of come aux party-to 'A lot of students came to the party' b. Indirect object: Dat.: Ugazabak langile asko-ri eskatu dio laguntza boss-erg worker a lot of-dat ask aux-sg help 'The boss asked a lot of workers to help' c. Direct object: Abs.: Mirenek liburu asko ikusi du liburutegian Miren-erg book a lot of see aux library-in 'Miren has seen a lot of books in the library'

  gazte anitzek altxatu behar izan zuen harria finally young many-erg lift must have aux-sg stone-D 'Ultimately, many youngsters had to lift the stone' Abs.: Ikasle anitz jin da festara student many come aux party-to b. Indirect object: Dat.: Ugazabak langile anitzi eskatu dio laguntza boss-erg worker many-dat ask aux-sg help 'The boss asked a lot of workers to help' c. Direct object: Abs.: Mirenek liburu anitz ikusi du liburutegian Miren-erg book many see aux library-in 'Miren has seen a lot of books in the library' Since non-agreeing quantifiers cannot take ergative case in this system (cf. (30a)), it is difficult to establish their distributive nature with the full set of tests we used in the previous cases. Take (30) and (31): (30) Azkenean gazte anitzek altxatu behar izan zuten harria finally young many-erg lift must have aux-pl stone-D 'Ultimately, many youngsters had to lift the stone' √ collective √ distributive (31) *Azkenean gazte anitzek altxatu behar izan zuen harria finally young many-erg lift must have aux-sg stone-D (33) a. Jonek liburu anitz ordenatu ditu Jon-erg book many arranged aux.pl 'Jon arranged many books' b. ??Jonek liburu anitz ordenatu du Jon-erg book many arranged aux.sg 'Jon arranged a lot of books' (34) a. Mikelek ikasle anitz ikusi ditu talde bakar bat osatzen Mikel-erg student many see aux.pl group single one forming 'Mikel has seen many students forming a single group' b. *Mikelek ikasle anitz ikusi du talde bakar bat osatzen Mikel-erg student many see aux.sg group single one forming 'Mikel has seen many students forming a single group'

(

  48) a. *Mendian ikusi ditut behi, ardi eta mando mountain-in see aux.pl cow sheep and mule b. Mendian ikusi ditut behiak, ardiak eta mandoak mountain-in see aux.pl cow-D.pl sheep-D.pl and mule-D.pl 'I have seen cows, sheeps, and mules in the mountain' (49) a. *Ricardok ardo edan du bazkaltzeko Ricardo-erg wine drink aux.sg lunch-for b. Ricardok ardoa edan du bazkaltzeko Ricardo-erg wine-D.sg drink aux.sg lunch-for 'Ricardo drank wine for lunch' In these dialects (i.e. Central-western and Transitional), the [DO+D] sequences in the examples in (48b) and (49b) can get two interpretations: (i) specific (equivalent to the reading we get by using French les, Spanish los or English the); (ii) existential (talking about a nonspecific set or quantity; see section 8.3) 8.3. Romance languages The distribution of Souletin BNs reminds us of the distribution of bare nouns in Romance languages (some of them, at least). In Spanish, for example, BNs are only possible in DO position (50a,b,c) -or in postverbal S position, (49d) -; cf. Bosque (1996) for a complete description of BNs in Spanish. (50) a. Juan ha comido patatas 'Juan has eaten potatoes' b. Pedro ha visto leones 'Pedro has seen lions' c. Mikel ha bebido café 'Mikel has drunk coffee' d. Llegaron estudiantes 'Students arrived' The BNs in (50) can only get the so-called existential interpretation, as was the case in Souletin. Now, in subjects of transitive predicates and subjects not allowing an existential reading, the presence of the determiner is necessary for the sentence to be grammatical. ix (51) a. *(Los) estudiantes han comido patatas the.pl students aux eaten potatoes 'The students has eaten potatoes' b. *(Los) dinosaurios están extintos the.pl dinosaurs are extinct 'Dinosaurs are extinct' c. *(Las) girafas son altas the.pl giraffes are tall 'Giraffes are tall' 8.4. Number marking: Spanish vs. Souletin Despite similarities, there's a significant difference between BNs in Souletin and in Romance languages: number marking. Spanish BNs must necessarily appear with the plural marker [-s] (52), whereas Souletin BNs are real BNs with no marker at all (53). absence of number marking in Souletin BNs does not eliminate their plurality (cf. (45)). x At this stage, a question that comes to our mind is the following: do Souletin BNs possess a non-overt plural marker? In Basque, the overt number marker is [-k].

  shows that masses cannot entertain a predication relation with non-homogeneous predicates. Homogeneous predicates are those that are both cumulative and divisive. The examples in (16) involve a nonhomogeneous predicate (to weigh more than 300 kilos). Whereas mass quantifications can not

	be the subject of the non-homogeneous predicate (16a), non-agreeing quantifiers with a count
	noun can (16b).
	(16) a. *Ur askok 300 kilo baino gehiago pisatzen du
	water a lot of 300 kilo than more weight-hab aux
	'*A lot of water weights more than 300 kilos'

b. Zaldi askok 300 kilo baino gehiago pisatzen du horse a lot of 300 kilo than more weight-hab aux 'A lot of horses weight more than 300 kilos'

  involves an agreeing vague quantifier. This yields two possible readings for the predicate: a distributive one, where each of the youngsters lifts the stone, and a collective one, where the entire set of youngsters lifts the stone. (19) also allows intermediate readings, where the set of youngsters divides in small groups to lift the stone. The range of distributive

	(26) *Ikasle bakoitzak elkar 'Jon has seen many students: Jon, Mikel, Pello, Martxel…' ikusi du
	(20) Azkenean gazte askok (24) a. Ikasle askok elkarren/bata bestearen altxatu behar izan zuen harria antz handia dute student each-erg reciprocal seen aux.sg b. *Jonek ikasle asko ikusi du: Jon, Mikel, Pello, Martxel…
	finally student many-erg each other-gen/one another-gen look big young a lot of-erg lift must have aux-sg stone-D '*Each student has seen each other' Jon-erg student a lot of see aux	aux.pl
	'Ultimately, many youngsters had to lift the stone' 'Many students look like each other/one another' 'Jon has seen a lot of students: Jon, Mikel, Pello, Martxel…'
	* collective b. *Ikasle askok elkarren/bata bestearen √ distributive 4.3. Enumeration and anaphora antz handia du
	student a lot of each other-gen/one another-gen look big	aux.sg
	'A lot of students look like each other/one another' (25) a. Irakasle askok elkar/bata bestea iraintzen dute 'Many professors insult each other/one another' b. *Irakasle askok elkar/bata bestea iraintzen du professor a lot of-erg each other/one another insult aux.sg 'A lot of professors insult each other/one another' Another difference between agreeing and non-agreeing quantifiers is that the latter cannot make reference to specific individuals. Thus, non-agreeing quantifiers cannot be antecedent to anaphoras, in opposition to what happens with agreeing quantifiers, as the examples in (27) (19) 4.2.2. Predicate classes professor many-erg each other/one another insult aux.pl show.
	We adopt Heim, Lasnik & May's (1991) analysis of reciprocals: in their view, reciprocals
	Non-agreeing quantifiers are incompatible with collective predicates (predicates that do not are complex quantificational expressions containing a distributive quantifier. This distributive
	allow event distribution). The examples in (21)-(23) all contain a predicate that does not quantifier is overt in some languages (cf. English each other). But if reciprocals possess a
	naturally allow atomic distribution (distribution down to the atomic entities making up a tacit distributive operator themselves, then the incompatibility between non-agreeing In (27a), we see that agreeing Basque quantifiers allow the enumeration of individuals, i.e.
	plurality). Whereas agreeing quantifiers can be combined with those predicates (a), non-quantifiers and reciprocals can be easily explained: the distributive operator requires a plural it is possible to make reference to the members of the set we are talking about. The
	agreeing ones cannot (b): set to operate on, one that can be broken into individual atoms. But if non-agreeing quantifiers enumeration of individuals denoted by the NP combined with non-agreeing quantifiers is not
	19) Azkenean gazte askok are themselves distributive, there is no plural set to operate on. The incompatibility between altxatu behar izan zuten harria possible, (27b).
	finally (21) a. Ikasle ohi askok festa horretan topo egin zuten young many-erg lifted must have aux-pl stone-D reciprocals and non-agreeing quantifiers is thus a subcase of 'vacuous quantification'. The
	'Ultimately, many youngsters had to lift the stone' student ex many-erg party that-in meet done aux-pl effect is analogous to (29), with a strong distributive quantifier: (28) a. Jonek ikasle asko ikusi ditu: Jon, Mikel, Pello, Martxel…
	√ collective 'Many ex-students met at that party' √ distributive Jon-erg student many see aux

readings in (

19

) is typical of count plural entities (see

Krifka, 1992)

. Unlike (

19

), (20) only allows a strict distributive reading, where youngsters individually lift the stone, and several stone-liftings (as many as there are youngsters) occur.

(27) a. Bezero asko i sartu dira gaur. _ i ez dira oso pozik atera. customer many come aux today neg aux very happy leave 'Many customers came today. They didn't leave very happy' b. *Bezero asko i sartu da gaur. _ i ez da oso pozik atera. customer a lot of come aux today neg aux very happy leave 'A lot of customers came today. They didn't leave very happy'

  prep), on the other hand, propose (in line with Borer 2005) that number is represented by two distinct syntactic positions in nominal expressions: by means of a number head ([-k] inBasque) and by means of a classifier head which portions-out count BNs to make them countable in order to interact with the counting function. In Central and Western dialects, this covert function is only available when a vague quantifier is present. In Souletin, it would seem that this covert function can apply directly on the noun. This classifier does not have phonological realization, but it is able to pluralize BNs in this dialect. Note that non-agreeing weak quantifiers show exactly the same syntactic distribution as BNs. We leave the relationship between these two phenomena for future research (cf.Etxeberria & Etxepare, in prep).

Acknowledgements

The research conducing to this paper has benefited from the Basque Government projects GIC07/144-IT-210-07 and Hm-2008-1-10, from the project FR2559 Fèderation Typologie et Universaux Linguistiques, from the project TSABL (ANR-07-CORP-033) from ANR, from the project FFI2008-00240 from MCE, from the project (MICINN)/FFI2011-29218 from MCE, from the project FFI2001-26906 from MCE as well as from the project ISQI from ANR. Thanks also to the two anonymous reviewers for useful comments and suggestions.

Usual disclaimers apply.

Reciprocals

Non-agreeing quantifiers, unlike agreeing ones, are incompatible with reciprocals:

'Ultimately, many youngsters had to lift the stone' * collective * distributive

The sentence in (30), with a agreeing weak quantifier in subject position, can obtain two interpretations, a collective one and a distributive one (just as was the case in the Centralwestern system). Now, the fact that non-agreeing quantifiers do not accept ergative case makes it impossible to conclude whether there are any differences in this respect between the Central-western and the Transitional system.

However, the next two subsections make it clear that non-agreeing quantifiers in the transition system are also distributive.

Predicate classes