Ugo Marchand
email: marchand@ircam.fr

Quentin Fresnel
email: fresnel@ircam.fr

Geoffroy Peeters
email: peeters@ircam.fr

GTZAN-RHYTHM: EXTENDING THE GTZAN TEST-SET WITH BEAT, DOWNBEAT AND SWING ANNOTATIONS

We present here the GTZAN-Rhythm test-set which is an extension of the GTZAN test-set with rhythm annotations. While there is a controversy related to the use of the GTZAN test-set for music genre classification [13], its audio content is however representative of real commercial music of various music genre and its audio is also already widely shared. Also, considering the various genres included in it, it provides a good balance between various types of rhythm including tracks with and without swing. We manually annotated the 1000 audio excerpts into beat, downbeat and swing annotations. We describe here this new test-set, as well as the validation of the annotations and the various uses one can make out of it (swing analysis, beat and downbeat tracking).

INTRODUCTION

The GTZAN-Rhythm test-set 1 is an extension of the GTZAN test-set with rhythm annotations (beat, downbeat and swing). While there exist several test-set related to rhythm annotations (see Table below), none of them provides annotation at the swing level (systematic displacement of the eight note).

We therefore decided to create a new test-set with the full metrical information including the swing factor for our research on swing estimation [START_REF] Marchand | Swing ratio estimation[END_REF]. Rather than starting from scratch, we have decided to extend an existing corpus, the GTZAN test-set, because its audio represents real commercial music and because its audio is already widely shared.

GTZAN-RHYTHM TEST-SET

In the following we describe this new test-set using the recommendation made by [START_REF] Peeters | Towards a (better) definition of the description of annotated m.i.r. corpora[END_REF] for "the description of annotated MIR corpora". The letters and numbers in brackets (such as "(B31)") refer to the description of [START_REF] Peeters | Towards a (better) definition of the description of annotated m.i.r. corpora[END_REF].

(T) - B - - ISMIR-2004 Songs yes T - - - - ISMIR-2004 Loop yes T - - - - ISMIR-2004 Ballroom [7] yes (T) (M) B D - Klapuri [6] on-demand (T) - B D - Hainsworth on-demand (T) - B - - Mazurka (MIREX) no (T) - B - - Isophonic-Beatles [10] yes (T) (M) B D - RWC [2] yes (T) (M) B D - SMC / MIREX-2012 [4] yes (T) - B - - Last-FM [8] yes T - - - - HJDB [3] no - - - D - GTZAN-Rhythm [9] yes T M B D S

(A) Raw corpus (audio):

The GTZAN-Rhythm test-set is based on the audio of a widely used test-set for music genre recognition research: the GTZAN test-set [START_REF] Tzanetakis | Musical genre classification of audio signals. Speech and Audio Processing[END_REF]. This test-set is made of 1000 audio excerpts of 30-seconds duration in mono at 22.050 Hz.

Those are equally distributed into 10 music genres: blues, classical, country, disco, hip-hop, jazz, metal, pop, reggae and rock.

(B) Annotations

(B1) Annotation origin: Each track of the test-set has been manually annotated (B15).

(B21) Concepts definition: The concepts being annotated are the downbeat, beat/tactus (quarter-note) and eight-note/tatum positions2 .

(B22) Annotation rules: The annotation was performed on the whole duration of each excerpt in a semiautomatic way. Beat and downbeat positions have been first automatically estimated using ircambeat software [START_REF] Peeters | Simultaneous beat and downbeat-tracking using a probabilistic framework: theory and largescale evaluation. Audio, Speech, and Language Processing[END_REF]. Then the resulting beat and downbeat positions have been manually corrected. The subdivision (by 2 or 3) of the beat positions have been used to generate the tatum positions which have been displaced by hand to fit locally the audio content.

(B31) Annotators: The annotations have been made by two different annotators, both are practicing musicians and researchers. One annotator did all the tracks numbered from 0 to 49, the other annotator did the tracks from 50 to 99. Each annotator then verified the work of each other, and reached a consensus when they were not in agreement.

(B32) Corpus Validation/Reliability: We tested the inter-annotator agreement between the two annotators. For this a part of the test-set (tracks from 95 to 99 from each genre) was annotated by both annotators.

For beat positions, we tested the inter-annotator agreement using the 3 followng metrics: F-measure, AMLt and Information Gain. We used the beat-tracking-evaluationtoolbox presented in [START_REF] Matthew Ep Davies | Evaluation methods for musical audio beat tracking algorithms[END_REF] to compute these metrics. The results are presented in Table 1. The good F-measure (91%), AMLt (95%) and information gain (4 bits) show that both annotators are most of the times in agreement for the positions of the beats. F-measure AMLt Inf. Gain Agreement 90.7% 94.9% 3.99 bits Table 1. Agreement between the annotators on beattracking, evaluated on a 5% subset of our base.

For downbeat positions, we tested the agreement for the subset of the most difficult tracks, the 153 tracks with swing. Knowing that the beat positons are shared by the two annotators, the two annotators agreed for the ones being the downbeat for 150/153 tracks. After discussion, they always reached a consensus.

(B4) Annotation tools: The annotation have been performed using the Audiosculpt 3.0 graphical software which allows to represent graphically the downbeats, beats and tatums positions over a spectrogram and allows dedicated editing of those marker types. (C2) Storage:

The GTZAN-Rhythm test-set can be downloaded at http://anasynth.ircam.fr/ home/media/GTZAN-rhythm/. The test-set contains:

• a set of 1000 musicdescription XML files (one for each audio track) containing the beat, downbeat and tatum positions,

• a single CSV file describing on each row the global information associated to a track: artist, title and a summary of the rhythm information (mean tempo, mean meter, mean swing ratio . . . ). The summary of the rhythm information has been obtained using the included python script generate.py. This script also exemplifies the reading of the musicdescription xml files,

• a set of 1000 JAMS3 files containing the beat, downbeat, tatum, artist, title and mean tempo, ... The artist and the title information have been kindly provided by Bob Sturm. If you use the GTZAN-Rhythm test-set, please cite [START_REF] Marchand | Swing ratio estimation[END_REF].

2. 3 (

 3 C) Documenting and storing (C1) Corpus Identifier: The identifier of the annotated corpus is GTZAN-Rhythm-V02.

Acknowledgements: This work was partly founded by the French government Programme Investissements d'Avenir (PIA) through the Bee Music Project.

We only annotated the eight-note position when swing exist or in case of a ternary track.

JSON Annotated Music Specification for Reproducible MIR Research[START_REF] Humphrey | Jams: a json annotated music specification for reproducible mir research[END_REF] 3. APPLICATIONSAs there is currently no test-set on which the swing ratio is annotated, our main contribution is to provide one. This test-set has been already used for our study on swing ratio estimation[START_REF] Marchand | Swing ratio estimation[END_REF]. We hope it could serve to quantitatively evaluate other related works on swing.

This test-set also provides beat and downbeat annotations on a large (1000 tracks) set of commercial audio tracks. It could therefore be used for the evaluation of beat-tracking, downbeat-tracking, tempo or meter estimation algorithms. Even though some audio tracks of the GTZAN-Rhythm test-set are quite simple for beat and downbeat tracking, some other (especially in the jazz, blues and classical section) are quite challenging.


R Triay

H H Fliche

Euler-Poisson-Newton approach in Cosmology

This lecture provides us with Newtonian approaches for the interpretation of two puzzling cosmological observations that are still discussed subject : a bulk flow and a foam like structure in the distribution of galaxies. For the first one, we model the motions describing all planar distortions from Hubble flow, in addition of two classes of planar-axial distortions with or without rotation, when spatial distribution of gravitational sources is homogenous. This provides us with an alternative to models which assume the presence of gravitational structures similar to Great Attractor as origin of a bulk flow. For the second one, the model accounts for an isotropic universe constituted by a spherical void surrounded by a uniform distribution of dust. It does not correspond to the usual embedding of a void solution into a cosmological background solution, but to a global solution of fluid mechanics. The general behavior of the void expansion shows a huge initial burst, which freezes asymptotically up to match Hubble expansion. While the corrective factor to Hubble law on the shell depends weakly on cosmological constant at early stages, it enables us to disentangle significantly cosmological models around redshift z ∼ 1.7. The magnification of spherical voids increases with the density parameter and with the cosmological constant. An interesting feature is that for spatially closed Friedmann models, the empty regions are swept out, what provides us with a stability criterion.

INTRODUCTION

On one hand, it is a matter of fact that the use of Newtons theory of gravity for understanding cosmology structures is easier than general relativity. On the other hand, it is however clear that such an approach must be taken with caution at large scales although it "is much closer to general relativity (GR) than commonly appreciated" [39]. With this in mind, we focus on Newtonian interpretation of two cosmological structures that have been mentioned in the literature, namely : the presence of a bulk flow and a foam like structure in the observed distribution of galaxies. For that, we write Euler-Poisson equations system in adapted coordinates, what enables us to define in a straightforward way solutions, as Newton-Friedmann and the vacuum (de Sitter) models. The bulk flow is analysed within an an-isotropic extension of Hubble model, see ([16]). Subsequently, the model of a single spherical void is derived within a covariant approach, see ( [17]).

EULER-POISSON EQUATIONS SYSTEM

As usual, for modelling the dynamics of the cosmological expansion, the distribution of gravitational sources is supposed to behave as dust (i.e., such a description does not account for shear, collision, etc. . . and the presence of radiations). The motions of such a pressureless medium are described at position r and time t by its specific density ρ = ρ( r,t) and velocity v = v( r,t). These fields are constraint by Euler equations system

∂ ρ ∂ t + div(ρ v) = 0 (1) ∂ v ∂ t + ∂ v ∂ r v = g (2)
where g = g( r,t) stands for the gravitational field. It satisfies the modified Poisson-Newton equations

-→ rot g = 0 (3) div g = -4πGρ + Λ ( 4
)
where G is Newton constant of gravitation and Λ the cosmological constant.

For the investigation of an expanding medium, it is convenient to write these equations with new coordinates, herein named reference coordinates. They are defined by

(t, x = r a ), a > 0 ( 5
)
where a = a(t) is a non constant monotonic function with a • = a(t • ) = 1. Such a method is motivated by the proper dynamics of medium, so that the equations system which accounts for the motion translates into its simplest form.

From now, let div and

-→
rot denote the differential operators with respect to x. Euler-Poisson-Newton (EPN) equations system reads in term of reference coordinates as follows

∂ ρ c ∂ t + div(ρ c v c ) = 0 (6) ∂ v c ∂ t + ∂ v c ∂ x v c + 2H v c = g c (7)
-→

rot g c = 0 (8)

div g c = - 4πG a 3 ρ c -3 Ḣ + H 2 - Λ 3 ( 9
)
where the dotted variables stand for time derivatives,

ρ c = ρa 3 , v c = d x dt = 1 a ( v -H r) , H = ȧ a ( 10
)
act respectively as the density and the velocity fields of medium in the reference frame, and

g c = 1 a g -Ḣ + H 2 r ( 11
)
as the acceleration field. Let us note that this model depends on two dimensionless parameter

Ω • = 8πG 3H 2 • ρ c , λ • = Λ 3H 2 • ( 12
)
where

H • = H(t • ).

Newton-Friedmann and Vacuum models

For these models, the most adapted family of functions a = a(t) verify a Friedmann type equation

Ḣ + H 2 - Λ 3 + 4πG 3 ρ • a -3 = 0 (13)
where ρ • is an arbitrary constant which stands for a parameter of function a(t). This differential equation admits the following integration constant

K • = 8πG 3 ρ • + Λ 3 -H 2 • ( 14
)
It integrates for providing us with the function H(t) from

H 2 = Λ 3 - K • a 2 + 8πG 3 ρ • a 3 ≥ 0 ( 15
)
Let us denote its asymptotical value

H ∞ = lim a→∞ H = Λ 3 ( 16
)
According to eq. ( 14), a depends on H • and on two additional parameters chosen among Λ, ρ • and K • . At this step, excepted Λ, these variables do not identify to cosmological parameters associated to Friedmann solution, they intervene merely in the coordinates choice. The formal expression t → a is derived as reciprocal mapping of a quadrature from eq. ( 15), it stands for the reference frame chronology.

It is obvious that EPN equation system eq. [START_REF] Conroy | [END_REF][START_REF] Courteau | Astrophysical Ages and Times Scales[END_REF][START_REF] Croton | The 2dFGRS team[END_REF]9) shows two trivial solutions :

1. The one defined by

ρ c = ρ • , v c = 0, g c = 0 ( 17
)
which accounts for a uniform distribution of dust, herein named Newton-Friedmann model (NF). According to eq. [START_REF] Fliche | Évaluation des paramètres cosmologiques à l'aide des propriétés optiques des quasars[END_REF]10), a and H recover now their usual interpretations in cosmology as expansion parameter and Hubble parameter respectively, ρ • = ρa 3 identifies to the density of sources in the comoving space and K • interprets in GR as its scalar curvature. We limit our investigation to motions which do not correspond to cosmological bouncing solutions, what requires the constraint

K 3 • < (4πGρ • ) 2 Λ (18)
to be fulfilled, according to analysis on roots of third degree polynomials. Hence, the kinematics shows two distinct behaviors characterized by the sign of K • . Indeed, H decreases with time by reaching H ∞ either upward (K • ≤ 0) or downward (K • > 0) from a minimum defined by

H m = H ∞ 1 - K 3 • Λ (4πGρ • ) 2 < H ∞ (19) at (epoch) a = 4πGρ • K -1
• , what defines a loitering period.

The one defined by

ρ c = 0, v c = (H ∞ -H) x, g c = 4πG a 3 ρ • x (20)
which accounts for vacuum solution, herein named Vacuum model (V).

ANISOTROPIC HUBBLE MODEL

This model accounts for collisionless motions which satisfy the following kinematics

r = A r • , A(t • ) = 1l ( 21
)
where A = A(t) stands for a non vanishing determinant matrix; its coefficients can be determined by an observer at rest with respect to Cosmological Background Radiation (CMB). One has detA > 0 because it identifies to a unit matrix at time t = t • . From eq. ( 21), the velocity field

v = v( r,t) reads v = d r dt = ȦA -1 r, Ȧ = d dt A (22)
where A -1 stands for the inverse matrix 3 . Let us choose

a(t) = 3 √ detA (23)
and with the unique matrix decomposition4

ȦA -1 = H1l + H • a 2 B, H = ȧ a , H • = H(t • ), tr B = 0 (24)
eq. (9,10,11) provide us with the equations of motion in the reference frame

v c = H • a 2 B x ( 25
)
g c = H • a 2 Ḃ + H • a 2 B 2 x ( 26
)
H 2 • a 4 trB 2 = - 4πG a 3 ρ c -3 Ḣ + H 2 - Λ 3 ( 27
)
where a stands for the (generalized) expansion factor, thus H = H(t) acts as the usual Hubble factor, and B = B(t) for a traceless matrix herein named distortion matrix. It characterises the deviation from isotropy of the (dimensionless) velocity field, its amplitude is defined by the matrix norm

B = tr (B t B) (28)
where the sign " t " stands for the matrix transposition. According to eq. ( 12), eq. ( 27) reads as follows

1 H 2 • Ḣ + H 2 = - Ω • 2a 3 + λ • - 1 3a 4 β 2 (29)
where

β n (t) = tr (B n ) , n = 1, 2, 3 (30)
Eq (29) shows that Ω • does not depend on spatial coordinates, i.e. the space distribution of gravitational sources is necessarily homogenous ρ = ρ(t) (31) and hence it stands for a constant parameter because of eq. [START_REF] Conroy | [END_REF]. By multiplying each term of eq. ( 29) by 2 ȧa, one identifies the following constant of motion

κ • = Ω • a + λ • a 2 - ȧ2 H 2 • - 2 3 a 1 β 2 a 3 da = Ω • + λ • -1 (32)
It is important to note that the chronology 5 , which is given by

dt = 1 H • ada P(a) ( 33
)
where

P(a) = λ • a 4 -κ • a 2 + Ω • a - 2 3 a 2 a 1 β 2 a 3 da ≥ 0, P(1) = 1 ( 34
)
is distortion dependent.

A Class of possible motions

Accordingly to usual Hubble model, if we assume a radial acceleration field

g c ∝ x ( 35
)
then the distortion matrix satisfies the evolution equation

Ḃ + H • a 2 B 2 = 1 3 H • a 2 trB 2 (36)
according to eq. [START_REF] Joevee | The Large Scale Structure of the Universe[END_REF][START_REF] Kirshner | [END_REF]. With

dτ = H • dt a 2 (37)
eq. ( 36) reads in a dimensionless form dB dτ

= 1 3 β 2 1l -B 2 (38)
With eq. (25,[START_REF] Joevee | The Large Scale Structure of the Universe[END_REF]33,36), the equations of motion read

d x dτ = B x, d 2 x dτ 2 = 1 3 β 2 x (39) dτ = da a P(a) ( 40
)
The resolution of these equations can be performed by mean of numerical techniques : the evolution with time of distortion matrix B is derived from eq. ( 38), the particles trajectories τ → x(τ) are obtained by integrating eq. ( 39) and the evolution of the generalized expansion factor a from eq. ( 40).

Analysis of analytic solutions

Analytic solutions of eq. (39,[START_REF] Saunders | Mapping the Hidden Universe: The Universe behind the Milky Way -The Universe in HI[END_REF]38) can be obtained thanks to particular properties of distortion matrix B. The parameters λ • , Ω • and H • given in eq. ( 12) correspond to cosmological parameters of Friedmann-Lemaître (FL) solution. The constraint β 2 = 0 in eq. (33,34) provides us with the FL chronology, where κ • in eq. ( 32) represents the curvature parameter in the FL model (i.e. the dimensionless scalar curvature Ω k of the comoving space, see [49]), while the flatness of (simultaneous events) Newton space. It must be noted that the particle position x as defined in eq. ( 5) does not identify to the usual FL comoving coordinate because the (generalized) expansion factor a depends on the anisotropy unless β 2 = 0.

Evolution of functions β n=2,3 . Because B is a traceless matrix, its characteristic polynomial reads

Q(s) = det (s1l -B) = s 3 - 1 2 β 2 s - 1 3 β 3 ( 41
)
according to Leverrier-Souriau's algorithm [START_REF] Souriau | Calcul linéaire[END_REF]. With Cayley-Hamilton's theorem (i.e. Q(B) = 0) and eq. ( 38) we obtain the following differential equations system

d dτ β 2 = -2β 3 ( 42
)
d dτ β 3 = - 1 2 β 2 2 ( 43
)
and we note that the discriminant of third order polynomial Q, it is proportional to

α = 3β 2 3 - 1 2 β 3 2 (44)
is a constant of motion (i.e., dα/dτ = 0). The integration of eq. (42,43) gives

τ = τ • + ε √ 6 2 β 2 (τ) β 2 (τ • ) dx √ 2α + x 3 , ε = ±1 (45)
Hence, β 2 is defined by a quadrature, and β 3 from eq. ( 42); in addition of the singular solution

β 2 = β 3 = 0, (i.e., B 3 = 0) (46)
defined equivalently either by β 2 = 0 or β 3 = 0, according to eq. (42,43).

The related dynamics depends on roots η i=1,2,3 of characteristic polynomial Q given in eq. ( 41), i.e. the eigenvalues of distortion matrix B. Their real values identify to dilatation rates at time τ toward the corresponding (time dependent) eigenvectors (not necessarily orthogonal). Their sum is null (β 1 = 0) and their product (β 3 = 3detB) is either decreasing with time or is null, according to eq. ( 43). The sign of α given in eq. ( 44) is used to classify the solutions as follows :

• if α = 0 then Q has a real double root η 1 = η 2 and a simple one η 3 . The related instantaneous kinematic shows a planar-axial symmetry (either a contraction within a plane with an expansion toward a transverse direction or vice versa), see sec. . If η 1 = η 3 then both vanish and the related solution identifies to the singular one defined in eq. ( 46); • if α > 0 then Q has a single real root η 1 ;

• if α < 0 then Q has three distinct real roots η i=1,2,3 . Their order is conserved during the evolution (since a coincidence of eigenvalues makes α = 0), the largest one must be positive while the smallest one must be negative (because β 1 = 0).

Planar kinematics. The singular solution B 3 = 0 shows a FL chronology and the distortion matrix

B = -B 2 • τ + B • , B 3 • = 0 (47)
is solely defined by its initial value B • , according to eq. ( 38). It is neither symmetric nor asymmetric (otherwise it vanishes), see eq. ( 28). Hence, eq. ( 39) transforms

d x dτ = -B 2 • τ + B • x ( 48
)
which accounts for eternal motions

x = exp -B 2 • τ 2 2 + B • τ x • = (1 + B • τ) x • (49)
The trajectory of a particle located at initial position x • identifies to a straight line toward the direction B • x • . The analysis of B • range (i.e., its image) provides us with characteristics of trajectories flow. The nilpotent property of B • shows that its kernel is not empty Ker(B • ) = / 0. Its dimension dim (Ker(B • )) = m characterizes the kinematics, which is either planar (m = 1) or directional (m = 2), i.e. a bulk flow. Conversely, if the kernel of distortion matrix B is not empty then β 3 = 3 det(B) = 0, and thus β 2 = 0, see to eq. (42,43). Therefore, all planar kinematics can be described by such a model.

Planar-Axial kinematics. If β 2 = 0 then the chronology differentiates from FL one. Let us focus on the α = 0 with two distinct eigenvalues η 1 = η 3 class of solutions. With eq. ( 45), eq. (42,43) integrate

β n = 6 (τ -τ ⋆ ) n , (n = 2, 3), τ ⋆ = τ • + √ 6β -1/2 2 (τ • ) (50)
which shows a singularity at date τ = τ ⋆ > 0 that splits the motion in two regimes τ < τ ⋆ and τ > τ ⋆ . The complete investigation of this singularity problem demands to solve an integro-differential equation, see eq. (33,34). The roots of

Q read η 1 = 1 (τ ⋆ -τ) , η 3 = -2η 1 ( 51
)
where η 1 stands for the double root. Among others, two class of solutions are defined by mean of a constant (time independent) matrix P, the projector associated to η 1 , see [START_REF] Souriau | Calcul linéaire[END_REF],

P 2 = P, trP = 2 (52)
They describe distinct kinematics depending on whether matrix B is diagonalizable.

• Irrotational motions : If B is diagonalizable then

B = η 1 (3P -21l) (53)
From eq. (39,53), one has

d x dτ = η 1 (3P -21l) x (54)
and the solution reads

x = -η 1 P ξ + 1 η 2 1 (1l -P) ξ ( 55
)
where ξ is constant. If the eigenvectors are orthogonal then the kinematics accounts for irrotational motions. • Rotational motions : If B is not diagonalizable then

B = η 1 (3P -21l) + 1 η 2 1 N (56) d x dτ = η 1 (3P -21l) + 1 η 2 1 N x ( 57
)
where N is a constant nilpotent matrix, which accounts for rotational motions on the eigenplane of P. The solution reads

x = -η 1 P ξ + 1 η 2 1 (1l -P) ξ + 1 η 2 1 N ξ ( 58
)
where ξ is constant.

Discussion

As a result, this anisotropic generalization requires an homogeneous distributions of matter. Because of the presence of strong density inhomogeneities in the sky distribution of galaxies catalogs, one is forced to ask whether it describes correctly the dynamics of observed cosmic structures. In principle, such a remark should be also sensible for questioning Hubble law when, regardless the isotropy, it is a fact that perturbations are not so dominant otherwise it would never have been highlighted. Actually, homogeneity is implicitly assumed in the standard cosmology for the interpretation of CMB isotropy and the redshift of distant sources, which provides us with an expanding background 6 . It is with such a schema in mind that this anisotropic Hubble law provides us with an hint on the behavior of the cosmic flow from decoupling era up to present date in order to answer whether the observed bulk flow is due exclusively to tidal forces.

Cosmic flow of flat LSS (B 3 = 0)

The B 3 = 0 class of solutions has interesting properties with regard to the stability of large scale structures that show a flat spatial distribution. To answer the question of whether observations define unambiguously the kinematics, the distortion matrix B is decomposed as follows

B = S + j( ω), trS = 0 ( 59
)
where S and j( ω) stand for its symmetric and its asymmetric7 component, and

ω = a 2 H • σ , σ = 1 2 -→ rot v (60)
accounts for the motion rotation, σ being the swirl vector. Hence, eq. ( 47) gives

B ω = S ω (61)
The evolution of the anisotropy with time is defined by

S = -S 2 • + j( ω • )j( ω • ) τ + S • (62) j( ω) = -(S • j( ω • ) + j( ω • )S • ) τ + j( ω • ) (63)
which couples the symmetric and the antisymmetric parts of the distortion matrix. The swirl magnitude reads

ω = ω, ω = 1 2 trS 2 (64)
according to eq. ( 60), since β 2 = 0. Its orientation cannot be determined from the data because the above equations describe two distinct kinematics corresponding to ± ω • that cannot be disentangle. According to eq. ( 59), if (and only if) the rotation ω = 0 then the distortion vanishes S = 0 since B is either a symmetric or antisymmetric. In other words, a planar distortion has necessarily to account for a rotation.

Constant distortion

Among above solutions which show planar kinematics, let us investigate the (simplest) one defined by B 2 • = 0. In such a case, k • ∝ ω and S ω = 0. According to eq. (62,63,64), linear calculus shows that the distortion is constant

S = S • , ω = ω • (65)
Such a distortion in the Hubble flows produces a rotating planar velocities field with magnitude ∝ H • a -2 . In the present case, the model parameters can be easily evaluated from data. The observed cosmic velocity fields are partially determined by their radial component

v r = v, r r = cz, r = r(m) = r u = a x, r = ct (66)
where m, z, u, t stand respectively for the apparent magnitude, the redshift, the line of sight, the photon emission date of the galaxy and c the speed of the light. According to eq. ( 22,24,59), the radial velocity of a galaxy located at position r is given by

v r = H + H u r, H u = H • a 2 u • S u (67)
Because trS = 0, it is clear that the sum of three radial velocities v r corresponding to galaxies located in the sky toward orthogonal directions and at same distance r provides us with the quantity H. Hence, simple algebra shows that the sample average of radial velocities within a sphere a radius r is equal to

v r = H r (68)
Therefore, for motions described by eq. ( 21), the statistics given in eq. ( 68) provides us with a genuine interpretation of Hubble parameter H = H(t). Hence, according to eq. (24,33,34), one obtains the (generalized) expansion factor

a(t) = exp t t • H(t)dt (69)
Hence, the cosmological parameters can be estimated by fitting the data to the function

ψ(t) λ • ,κ • ,Ω • = P(a) = a 2 H/H • , λ • + κ • + Ω • = 1 (70)
The component of matrix S can be estimated by substituting H in eq. ( 67), and ω is obtained from eq. ( 65). It is clear that the above model is derived in the Galilean reference frame, where the Euler-Poisson equations system can be applied. Hence, a non vanishing velocity of the observer with respect to this frame an produces a bipolar harmonic signal in the H u distribution of data in the sky, which can be (identified and then) subtracted.

Discussion

The dynamics of a homogenous medium and anisotropic moving under Newton gravity was already studied by describing the evolution of an ellipsoid [START_REF] Ya | The Structure and Evolution of the Universe[END_REF]. The current approach enables us to identify characteristics of the dynamics of the deformation from isotropic Hubble law in a more systematic way by mean of the distortion matrix.

At first glance, if the planar anisotropic of the space distribution of galaxies within the Local Super Cluster (LSC) is stable then the above solution can be used for understanding its cosmic velocity fields, k • being orthogonal to LSC plane. It is well known however that the distribution of galaxies is not so homogenous as that, whereas this model describes motions of an homogenous distribution of gravitational sources. On the other hand, such an approximation level is similar to the one which provides us with the observed Hubble law, that is included in this model (S = 0).

SPHERICAL VOIDS IN NEWTON-FRIEDMAN UNIVERSE

For modelling the dynamics of a spherical void in a uniform dust distribution, we use a covariant formulation of Euler-Poisson equations system (see Appendix). The model is obtained by sticking together the local solutions V and NF of Euler-Poisson equations system as defined previously in the Newton-Friedmann and Vacuum models, where the function a stands for the Friedmann expansion parameter. Such a model accounts for the dynamics of their common border (i.e. boundaries conditions), which is a material shell. For convenience in writing, the symbols S, V and NF denotes both the medium and the related dynamical model. A qualitative analysis of solutions is performed and a general discussion is given subsequently

Dynamical model

We consider three distinct media : a material shell with null thickness (S), an empty inside (V) and outside a uniform dust distribution (NF). These media behave such that S makes the juncture of V with NF as given by eq. (17,20). The tension-stress on S is assumed to be negligible, what is characterized by a (symmetric contravariant) mass-momentum tensor defined as follows

T 00 S = (ρ S ) c , T 0 j S = (ρ S ) c v j c , T jk S = (ρ S ) c v j c v k c ( 71
)
The background is described by the following mass-momentum tensor

T 00 NF = ρ c , T 0 j NF = 0, T jk NF = 0 (72)
According to Appendix , since the eulerian function(al)

T (x → γ) = T µν S γ µν dtdS + T µν NF γ µν dtdV ( 73
)
vanishes when γ reads in the form γ µν = 1 2 ∂µ ξ ν + ∂ν ξ µ , one has

S ∂ 0 ξ 0 + g j c ξ j + (∂ j ξ 0 + ∂ 0 ξ j -2Hξ j ) v j c + v j c v k c ∂ j ξ k (ρ S ) c dtx 2 dΩ = - NF ρ c ∂ 0 ξ 0 dtx 2 dxdΩ ( 74
)
where dΩ stands for the solid angle element. The radial symmetry of solutions enables us to write the reduced peculiar velocity and acceleration of a test particle located on the shell as follows

v c = α x, g c = β x (75)
where the functions α = α(t) and β = β (t) have to be determined. A by part integration of eq. ( 74) provides us with

t 2 t 1 (∂ 0 (ρ S ) c + 3(ρ S ) c α -ρ c αx) x 2 ξ 0 dt (76) = t 2 t 1 ∂ 0 ((ρ S ) c α) + 4(ρ S ) c α 2 + 2H(ρ S ) c α -(ρ S ) c β x -1 x 3 ξ dt
where x = x stands for the radius of S and ξ = ξ 2 1 + ξ 2 2 + ξ 2 2 . This equality must be fulfilled for all bounded time interval and compact support 1-form. Hence, we easily derive the conservation equations for the mass

∂ 0 (ρ S ) c + (3(ρ S ) c -ρ c x) α = 0 (77)
and for the momentum dα dt

+ 1 + ρ c (ρ S ) c x α 2 + 2Hα + β x = 0 ( 78
)
With eq. ( 5), the calculation of the gravitational force from the entire shell acting on a particular point 8 provides us with

β = 4πG a 3 ρ c 3 - (ρ S ) c 2x ( 79
)
About mass conservation, it is noticeable that

(ρ S ) c = 1 3 ρ c x (80)
is solution of eq. ( 77), what ensures that the amount of matter which forms the shell comes from its interior. Hence, eq. ( 78) transforms dα dt

+ 4α 2 + 2Hα - 2πG 3 ρ c a 3 = 0 (81)
It is convenient to use the dimensionless variable

χ = 4 α H • a 2 (82)
where the ratio α H -1

• stands for the expansion rate of S in the reference frame. Hence, eq. ( 81) transforms into a Riccati equation dχ da

= Ω • - χ 2 a 1 P(a) ( 83
)
where

P(a) = λ • a 4 -k • a 2 + Ω • a, k • = K • H 2 • , P(1) = 1 ( 84
)
where the dimensionless parameters 9 are defined in eq. [START_REF] Pj | Proc. of Ninth Marcel Grossmann Meeting on General Relativity[END_REF]14). According to eq. (10,82), the evolution of the radius of S is given by

x = x i exp a a i χda 4a P(a) ( 85
)
where x i and a i stands for the initial values at time t i . 8 The modified newtonian gravitation field reads g = Λ 3 -Gm r 3 r 9 These notations are preferred to the usual Ω Λ = λ • and Ω K = -k • for avoiding ambiguities on the interpretation of cosmological parameters, see e.g., [START_REF] Fliche | [END_REF], [14].

Qualitative analysis

The shell expansion is analysed in term of dimensionless quantities : the magnification X and the expansion rate Y

X = x x i , Y = α H • (86)
Their evolution with the expansion parameter a are obtained from eq. (83, 85) by numerical integration 10 . Let us focus our investigation around the generally accepted values λ • = 0.7 and Ω • = 0.3. The initial conditions lie on the expansion rate Y i and the formation date t i , as expressed by means of a i = a(t i ). The values a i = 0.003 and Y i = 0 (void initially expanding with Hubble flow) are used as standard in our discussion. Hereafter, we simply provides us with synthesis results, more information can be found in [17].

Kinematics

The expansion velocity of the shell S with respect to its centre reads

v = yH r, y = 1 + Y h , h = H H • (87)
where y stands for the corrective factor to Hubble expansion. The general trend of the kinematics can be derived from the diagram y versus a in Fig. 1. It results from the initial conditions a i = 0.003, Y i = 0, and the cosmological parameters

Ω • = 0.3. and λ • = 0, λ • = 0.7, λ • = 1.4.
The shell S expands faster than Hubble expansion (y > 1) at early stages of its evolution with no significant dependence on Λ. Λ-effect appears later by preserving the curve from an earlier decreasing. If λ • > 0.7 (i.e., k • > 0) then it is characterized by a significant protuberance at z ∼ 1.7, the larger the λ • (i.e., k • ) the higher the bump. It is due to the existence of a minimum value H m of Hubble parameter H which is reached during the cosmological expansion, see eq. (19). After this period the expansion of S reaches asymptotically Hubble behavior. It is interesting to note that the present epoch (a = 1) appears quite peculiar because of the relative proximity of curves, but it is solely an artifact 11 .

Dichotomy between cosmological parameters

With the aim of disentangling the effect of the cosmological constant on the evolution of spherical voids from that of the outer density (herein named Λ-effect and Ω-effect respectively), the dependence of X and Y on cosmological parameters Ω • and λ • are investigated separately within acceptable values domains in accordance with observations. 1. Ω-effect. -With a constant λ • = 0.7, the magnification X and the expansion rate Y increase with Ω • , what interprets by the attraction of shell particles toward denser regions. The growth shows a huge burst which freezes asymptotically up to matching Hubble expansion (Y = 0). This trend is not significantly modified by other acceptable values of λ • . 2. Λ-effect . -With a constant Ω • = 0.3, the larger the λ • the larger the magnification X. As for Ω-effect, the magnification X increases nonlinearly with λ • . With Λ > 0, this phenomenon interprets as a repulsive effect of gravity in empty regions (also named gravitational repulsion of vacuum). The expansion rate Y is weakly dependent on λ • , what means that it does not characterizes Λ. However, it has a cumulative late effect which is reflected on X. 3. Interpretation of parameter k • . -According to eq. ( 20), the test-particle within S moves in the reference frame toward its centre but if k • > 0 then it starts moving toward its border at date a(t ⋆ ) = 8 3 πGρ • K -1 • . Such a property of sweeping out the void region (due to H < H ∞ ) interprets as a stability criterion for void regions. A dimensional analysis of eq. ( 14) shows that the Newtonian interpretation of k • = λ • + Ω • -1 corresponds to a dimensionless binding energy for the universe 12 . Because K • is an integration constant, this provides us with a meaningful procedure for comparing the magnitudes of Λ-effect or Ω-effect between different world models at constant k • . For Friedmann model, k • stands for the dimensionless spatial curvature (of the comoving space), k • > 0 world models correspond to spatially closed universes. 4. Dependence on initial conditions. -The earlier the date of birth the larger the magnification X. One has a limiting curve a → Y defined by a i = 0, which characterizes Ω • . All expansion rate curves related to other formation dates are located in lower part and reach it asymptotically. One can expect that the effects on the dynamics of growth resulting from reasonable initial expansion rates Y i = 0, related to physical processes (e.g., supernovae explosions), are negligible at primordial epoch 13 , what legitimizes the initial condition Y i = 0 at a i = 0.003. Evolutions with other initial conditions on Y i can be deduced from that since any given point (a i ,Y i ) in the diagram belongs to a single evolution curve.

Synthesis

As a result, the expansion of a spherical void does not show a linear regime, but a huge initial burst which freezes asymptotically up to matching Hubble expansion. Its radius (in the reference frame) x increases from an initial size x i with the cosmological parameters Ω • and λ • . The related individual effects interpret respectively by the gravitational attraction from the outer parts and by the gravitational repulsion of its borders (or the vacuum from the inner parts). The larger these parameters the higher the magnification X = x/x i . The dynamics is sensitive to Ω-effect at primordial epochs and to Λ-effect later on. The evolution of its expansion velocity (in the reference frame) v c = Y H • x with time does not characterise λ • but Ω • . On the other hand, the cosmological constant intervenes significantly on the kinematics by means of the corrective factor y to Hubble expansion by preserving the expansion from an earlier decreasing. Moreover, a values domain of cosmological parameters, which corresponds to relativistic spatially closed world models, are characterised by a significant expansion of voids. The later reaches a maximum at redshift z ∼ 1.7 (with Ω • ∼ 0.3), the larger the λ • (i.e., k • ) the higher the expansion rate. Notwithstanding, the perturbation on redshift of sources located on the shell of expanding voids does not exceed

∆ z = XY 1 + z x i H • c ( 88
)
which is a tiny value (∼ 10 -3 ) because of counterbalancing behaviors of X and Y . The interior of voids shows a de Sitter expansion v = Λ/3 r which sweeps out it if k • > 0, what interprets as a stability criterion. The expansion rate evolutions related to formation dates show a common envelope curve which characterizes the density parameter Ω • .

Discussion

It is clear that such a simple model should be applied with caution to observed distribution of galaxies since the effects related to foamlike patterns are not considered. Moreover, because Newtonian approach is used instead of general relativity, it applies solely to reasonable sized voids 14 so that the main features should remain qualitatively reliable. On the other hand, based on the connection of two exact solutions of Euler-Poisson equations system, which is a huge advantage for investigating properly the behavior of a single void, it provides us with an appreciable hint to the dynamics with respect to cosmological parameters.

CONCLUSION

Euler-Poisson equations system written in an adapted coordinates system enables us to define in a straightforward way cosmological solutions, such as Newton-Friedmann and the vacuum (de Sitter) models, to investigate an anisotropic extension of Hubble flow and the dynamics of spherical voids in an expanding universe, as long as Newton approximation with cosmological constant is valid. The related results are summarised as follows :

• An-isotropic Hubble flow. The present solution generalizes the Hubble law and provides us with a better understanding of cosmic velocity fields within large scale structures. As a result, it implies necessarily an homogenous distribution of gravitational sources, as similarly to Hubble law. Because the chronology identifies to FL chronology for a vanishing distortion, this model interprets as a Newton approximation of anisotropic cosmological solutions. The motions are characterised by means of a constant of motion α. Among them, particular solutions can be easily derived for α = 0. They describe all planar distortions, in addition of two classes of planar-axial distortions with or without rotation. Among these solutions, the one which ensures a planar kinematics is of particular interest because it describes constant (eternal) and rotational distortions. This solution can be fully determined from observational data except for the orientation of the rotation. The sensible result is that the velocity field is not potential. It is interesting to note that this model accounts for motions which might be interpreted as due to tidal forces whereas the density is homogeneous. It is an alternative to models which assume the presence of gravitational structures similar to Great Attractor as origin of a bulk flow. • Spherical voids. The model of an isotropic universe constituted by a spherical void surrounded by a uniform distribution of dust is defined such that the behaviors of bordering regions inside and outside the shell, coincide with Newton-Friedmann (NF) solutions. The connection conditions between these two regions are investigated in Classical Mechanics, what provides us with the dynamics of the shell. Such a schema does not correspond to the usual embedding of a void solution into a cosmological background solution, but interprets as a non linear perturbation of NF models. The general behavior of the void expansion shows a huge initial burst, which freezes asymptotically up to match Hubble expansion. While the corrective factor to Hubble law on the shell depends weakly on cosmological constant at early stages, it enables us to disentangle significantly cosmological models around redshift z ∼ 1.7. The magnification of spherical voids increases with the density parameter and with the cosmological constant. An interesting feature is that for spatially closed Friedmann models, the empty regions are swept out, what provides us with a stability criterion.

APPENDIX : COVARIANT FORMULATION OF EULER EQUATIONS SYSTEM

Souriau's covariant formulation of Euler-Poisson equations [44] can be summarized as follows : The geometrical interpretation of Newton dynamics [11] shows that the component of gravitational field g identify to the only non 14 At first glance, on might argue that the related scale should be larger than the size of background structures in order to make the uniform density hypothesis an acceptable approximation, what looks like a dilemma. However, such an hypothesis, which applies to the entire universe, is ensured by the isotropy of CMB.

null Christoffel components of newtonian connection Γ j 00 = -g j Hence, one obtains their expression in the new coordinates system defined in eq. ( 5), and the only non null components read (Γ c ) where ∂ stands for the covariant derivative and x → ξ for a compact support 1-form. In such a case, it is obvious to show that ∂µ T µν = 0 interprets as Euler equations, in any coordinates systems.

FIGURE 1 .

 1 FIGURE 1. The corrective factor y to Hubble expansion.

c 4 T 4 T

 44 Let T be a function(al) defined on the set of symmetric density tensors x → γ on newtonian spacetime R (x → γ) = R µν γ µν dtdV where dtdV stands for the volume element and T for a symmetric contravariant tensors which accounts for the media. The measure density T is eulerian if and only if it vanishes for all covariant tensor fields which reads γ µν = 1 2 ∂µ ξ ν + ∂ν ξ µ

Unité Mixte de Recherche (UMR 6207) du CNRS, et des universités Aix-Marseille I, Aix-Marseille II et du Sud Toulon-Var. Laboratoire affilié à la FRUMAM (FR

2291).

A -1 A = AA -1 = 1l

By using the trace, the determinant one has d dt detA = Tr ȦA -1 detA.

According to observations, one assumes H • > 0 since the case H • < 0 accounts for a collapse, what is not envisaged.

Namely the comoving space of FL world model onto which the gravitational instability theory is applied for understanding the formation of cosmic structures.

The operator j stands for the vector product, u × ω = j( u)( ω).

Because the mapping t → a(t) = 1/(1 + z) is a monotonic function in the present investigation, their evolution with cosmic time t or with redshift z can be straightforwardly derived.

The lower k • the faster the cosmological expansion. Note that it works with an opposite direction for the expansion of spherical voids.

Since Hubble expansion is all the more important towards the past, the earlier the formation date the weaker this effect, according to eq. (87).

ACKNOWLEDGMENTS

The authors thanks warmly Mário Novello, and co-organizers, for this high level School, the hospitality and financial supports.


Sergio Peignier

Christophe Rigotti

Guillaume Beslon

Subspace Clustering for all Seasons

Subspace clustering is recognized as a more general and difficult task than standard clustering since it requires to identify not only objects sharing similar feature values but also the various subspaces where these similarities appear. Many approaches have been investigated for subspace clustering in the literature using various clustering paradigms. The reader is referred for instance to [START_REF] Müller | Evaluating clustering in subspace projections of high dimensional data[END_REF] for detailed reviews and comparisons of the best methods and main categories. Even if many evolutionary clustering approaches exist [START_REF] Hruschka | A survey of evolutionary algorithms for clustering[END_REF] very few of them address the subspace clustering problem and still include non-evolutionary stages [START_REF] Sarafis | Towards effective subspace clustering with an evolutionary algorithm[END_REF][START_REF] Vahdat | Bottom-up evolutionary subspace clustering[END_REF].

According to [START_REF] Banzhaf | Guidelines: From artificial evolution to computational evolution: a research agenda[END_REF], bio-inspired optimization algorithms could be improved by incorporating knowledge from molecular and evolutionary biology. A promising source of advances in optimization is one of the important phenomena in evolutionary biology: the dynamic evolution of the genome structure. Several studies showed for instance that an evolvable genome structure allows evolution to modify the effects that evolution operators (e.g., mutations) have on individuals, a phenomenon known as evolution of evolution [START_REF] Hindré | New insights into bacterial adaptation through in vivo and in silico experimental evolution[END_REF]. In this paper, we present Chameleoclust, an evolutionary subspace clustering algorithm that incorporates a genome having an evolvable structure. The genome is a coarsegrained genome, inspired on [START_REF] Crombach | Chromosome rearrangements and the evolution of genome structuring and adaptability[END_REF], and defined as a list of tuples (the "genes"), each tuple containing numbers. These tuples are mapped at the phenotype level to denote core point locations in different dimensions, which are then used to collectively build the subspace clusters, by grouping the data around the core points. The biological analogy here would be that each gene codes for a molecular product and that the combination of molecular products associated together codes for a function, i.e., a cluster. To allow for evolution of evolution, Chameleoclust genome contains a variable proportion of functional elements as , and is subject to local mutations and to large random rearrangements, namely: large deletions, duplications and translocations. Local mutations and rearrangements may thus modify the genome elements but also the genome structure. The key intuition in the design of the Chameleoclust algorithm is to take advantage of such an evolvable structure to detect various numbers of clusters in subspaces of various dimensions and to self-tune the main evolutionary parameters (e.g., levels of variability). Figure 1 illustrates Chameleoclust genome structure (genome length and proportion of functional elements) and fitness convergence on a synthetic benchmark dataset. The reader is refered to [START_REF] Peignier | Subspace clustering using evolvable genome structure[END_REF] for a detailed description of Chameleoclust. The algorithm has been assessed using a reference framework for subspace clustering evaluation, and compared to state-of-the-art algorithms on both real and synthetic datasets [START_REF] Müller | Evaluating clustering in subspace projections of high dimensional data[END_REF]. However in [START_REF] Müller | Evaluating clustering in subspace projections of high dimensional data[END_REF] the parameters of the different state-of-the-art algorithms have been optimized according to external in order to compare the different algorithms by the best possible subspace clusterings they could achieve if we were able to find the most appropriated parameter values. Since generally no external labeling is available when we search for clusters, parameter tuning is most of the time a difficult task and these high quality subspace clusterings are likely to be hard to obtain. This is one of the benefits evolution of evolution may offer. Consequently we decided not to perform any parameter optimization using external information for Chameleoclust. So, we compare clusterings effectively found by Chameleoclust to the best clusterings that could potentially be found by the other algorithms. every dataset, the performances of Chameleoclust are reasonable with respect to the best possible runs of the other algorithms (results on the real world dataset "Shape" are presented as an example in Figure 2). For all datasets, the number of clusters and the subspaces dimensionalities found by Chameleoclust are coherent with the number of clusters found by the other algorithms (e.g., Figure 2). Indeed Chameleoclust is able to adapt the number of clusters it produces and also their average dimensionality according to each dataset without specific parameter tuning, as summarized in Figure 2. In addition, we give in Figure 3 the runtime of Chameleoclust with respect to the number of objects and number of dimensions of the synthetic datasets. The corresponding curves show that the algorithm scales rather linearly in both cases. q q q q q q q 2000 3000 4000 5000 0 500 1000 2000 DATABASE SIZE RUNTIME (SECONDS) q q q q q Figure 3: Runtime. Maximal in red (circles), median in cyan (triangles), average in green (squares) and minimal in blue (crosses).

The results obtained with the Chameleoclust algorithm show that evolution of evolution, through an evolvable genome structure, is usefull to solve a difficult problem such as subspace clustering.

Figure 1 :

 1 Figure1: Average of the best fitness, genome length and proportion of functional elements (functionality ratio) of all individuals for 10 runs on a synthetic dataset (1500 objects, 10 dimensions, 0% noise).

Figure 2 :

 2 Figure 2: Average number of clusters and average dimensionality per cluster found for each real world dataset (left). Results for the Shape real dataset: 17 dimensions, 9 classes, 160 objects (right)

  For almost

	Dataset N umClusters AvgDim		Accuracy	CE	N umClusters AvgDim
	breast	15.2	5.32	CLIQUE	max min max min max 0.76 0.76 0.01 0.01 486	min 486	max min 3.3 3.3
	diabetes glass	60.8 34.8	2.06 2.67	DOC MINECLUS 0.79 0.60 0.58 0.46 64 0.79 0.54 0.56 0.38 53 SCHISM 0.74 0.49 0.10 0.00 8835	29 32 90	13.8 12.8 17.0 17.0 6.0 3.9
	liver	62.8	1.93	SUBCLU FIRES	0.70 0.64 0.00 0.00 3468 3337 0.51 0.44 0.20 0.13 10 5	4.5 4.1 7.6 5.3
	pendigits shape	17.5 14.9	5.61 7.27	INSCY PROCLUS P3C	0.76 0.48 0.18 0.16 185 0.72 0.71 0.25 0.18 34 0.61 0.61 0.14 0.14 9	48 34 9	9.8 9.5 13.0 7.0 4.1 4.1
	vowel	48.7	2.59	STATPC Chameleoclust 0.82 0.73 0.42 0.34 16 0.74 0.74 0.45 0.45 9	9 13	17.0 17.0 10 7.29

Acknowledgments This research has been supported by EU-FET grant EvoEvo (ICT-610427).


Charles Rocabert
email: charles.rocabert@inria.fr

Carole Knibbe

Guillaume Beslon

Towards a Integrated Evolutionary Model to Study Evolution of Evolution

Keywords: Evolution of Evolution, Robustness, Evolvability, Genotypeto-Phenotype Mapping, In silico Experimental Evolution

Variation and Selection are the two core processes of Darwinian Evolution. Yet, both are directly regulated by many processes that are themselves products of evolution. Microorganisms efficiently exploit this ability to dynamically adapt to new conditions. Thus, evolution seems to have optimised their own ability to evolve, as a primary mean to react to environmental changes. We call this process Evolution of Evolution (EvoEvo). In this paper, we propose to use an integrated evolutionary model including complex and evolvable genotype-to-phenotype mapping to study EvoEvo. As exemplified in this paper, the integrated evolutionary model will allow us to decipher the EvoEvo strategies and to offer new hypothesis and predictions on the evolution of microorganisms.

Introduction

Life on Earth evolved for billion years in ever changing environments, undergoing smooth or brutal, cyclic or unseen variations. To survive in such conditions, being adapted to the current environment is not enough. Extant organisms had to deal with the evolutionary competition but they also had to deal with the variations of their environment to stay adapted despite the rapid and sometimes profound crisis they had to cope with. How did they do so is an open question. Did the extant organisms survive by chance or did they survive because, being regularly confronted to such crisis, they evolved reaction/adaptation mechanisms?

Experimental evolution, where fast replicating organisms (e.g. bacteria or viruses) are evolved in controlled environments for thousands of generations [START_REF] Philippe | Evolution of global regulatory networks during a long term experiment with Escherichia coli[END_REF], allows to recover precisely the evolution history of lab strains by reviving frozen samples and performing data analysis. Many evolution results have shown that microorganisms are able to evolve at an amazing speed: in virtually all experimental frameworks that use bacteria or viruses, important phenotypic innovations have emerged in only a few tens of generations [START_REF] Rainey | Adaptive radiation in a heterogeneous environment[END_REF][START_REF] Zhang | Acceleration of emergence of bacterial antibiotic resistance in connected microenvironments[END_REF] and in many cases, evolution tend to be partly reproducible [START_REF] Hekstra | Contingency and statistical laws in replicate microbial closed ecosystems[END_REF][START_REF] Raeside | Large chromosomal rearrangements during a long-term evolution experiment with Escherichia coli[END_REF]. Microorganisms efficiently use mutation and selection to dynamically adapt to new conditions. Thus, evolution seems to have optimized their own ability to evolve, as a primary mean to react to environmental changes. We call this process "Evolution of Evolution" (EvoEvo, [START_REF] Hindr | New insights into bacterial adaptation through in vivo and in silico experimental evolution[END_REF]).

Experimental evolution, despite its explanatory and statistical power, remains a long and costly process. An alternative is to simulate evolution in a computer. However, Evolution of Evolution implies the interaction of a wide range of biological structures (e.g. genome, genetic regulation network, metabolic network, ...), so we need to develop complex models ("The world is complex and we need all the tools that we can muster to understand it" [START_REF] Peck | Simulation as experiment: a philosophical reassessment for biological modeling[END_REF]). Following this idea, in silico experimental evolution is a growing field in evolutionary biology (reviewed in [START_REF] Mozhayskiy | Microbial evolution in vivo and in silico: methods and applications[END_REF] and [START_REF] Hindr | New insights into bacterial adaptation through in vivo and in silico experimental evolution[END_REF]), and a lot of theoretical questions have been deciphered with this approach (e.g. [START_REF] Wilke | Evolution of digital organisms at high mutation rates leads to survival of the flattest[END_REF][START_REF] Knibbe | A long-term evolutionary pressure on the amount of noncoding DNA[END_REF][START_REF] Beslon | From digital genetics to knowledge discovery: Perspectives in genetic network understanding[END_REF][START_REF] Crombach | Evolution of evolvability in gene regulatory networks[END_REF][START_REF] Cuypers | Virtual genomes in flux: an interplay of neutrality and adaptability explains genome expansion and streamlining[END_REF]).

In this paper, we present a new tool to study Evolution of Evolution. we first remind the theoretical background of EvoEvo and why we should use an integrated evolutionary model to study it. In a second part, we will present our model, and finally, our working plan and the perspectives of the integrated evolutionary model.

What is Evolution of Evolution ?

Variation and Selection are the two core engines of Darwinian Evolution. Yet, both are directly regulated by many processes that are themselves products of evolution (e.g. DNA repair, mutator genes, transposable elements, horizontal transfer, stochasticity of gene expression, sex, network modularity, niche construction...). This ability is the core of the EvoEvo process. Different "Evolution of Evolution" strategies have been proposed in the literature, including regulation of variability, robustness/evolvability strategies, bet-hedging ... However, most of them are poorly characterized and the conditions under which they evolve as well as their consequences are generally unknown.

To develop a better understanding of EvoEvo, we propose to build a conceptual framework based on two important concepts of evolutionary biology: the genotype-to-phenotype mapping and the fitness landscape. The genotypeto-phenotype mapping summarizes in a single conceptual entity the complex molecular processes by which information flows from the genetic sequence to the organism's phenotype. It thus concatenates in a single abstract process different phenomena such as mRNA transcription, gene translation, protein folding, biochemistry and cell dynamics. The central idea of the fitness landscape is that organisms or populations in evolution can be represented as points on a landscape where the altitude represents the fitness. Selection can be represented by the local gradient of altitude and the mutation can be represented as a random noise added to individual positions.

We propose to study EvoEvo by focusing on four characteristics of the genotype-to-phenotype mapping and the fitness landscape:

Variability. Variability is the ability to generate new phenotypes, by mutations or by stochastic fluctuations. It is a necessary condition for any evolu-tionary process to take place. However, in biological organisms, the amount of variability is controlled by complex pathways that e.g. correct DNA mismatches or double-strand breaks. Mutational operators are highly diversified, including point mutations, but also large chromosomal rearrangements that can rapidly reshuffle the chromosome organisation, extend or reduce the gene repertoire of an organism or even duplicate its entire genome through whole genome duplication.

Robustness. Although mandatory, variability is a very dangerous process since it permanently produces deleterious mutations that lead to poorly adapted individuals. Robustness may evolve to correct these deleterious effects. It enables evolving systems to support mutational events without loosing fitness through e.g. canalisation or the selection of structures that creates neutral landscapes.

Evolvability. Depending on the genotype-to-phenotype mapping, the proportion of deleterious/neutral/favorable mutational events may change. Evolvability is the ability of a specific genotype-to-phenotype mapping to increase the proportion of favorable events. This can be done by the selection of specific genome structures or by the selection of specific network structures.

Open-endedness. Biological evolution is not directed towards a specific target. On the opposite, evolution has the ability to generate new challenges while evolving by e.g. exploiting new niches created by the evolution of other species.

The central concept of EvoEvo is the following: if the genotype-to-phenotype mapping and the fitness landscape are allowed to change over time, if they can be (indirectly) selected, then they can evolve and acquire properties than could favour evolution in changing environments.

Why an integrated evolutionary model ?

Computational models have been used to study evolution since the beginning of the 90th [START_REF] Adami | Digital genetics: unravelling the genetic basis of evolution[END_REF]. However, since then, most computational models used a partial representation of the genotype-to-phenotype mapping, generally in a fixed, predefined, fitness landscape. By simulating the evolution of such and such organisation level (the genome, the genetic regulation network, the metabolic network, ...), different authors have studied evolution of robustness, evolvability or variability of these specific levels [START_REF] Wilke | Evolution of digital organisms at high mutation rates leads to survival of the flattest[END_REF][START_REF] Knibbe | A long-term evolutionary pressure on the amount of noncoding DNA[END_REF][START_REF] Crombach | Evolution of evolvability in gene regulatory networks[END_REF][START_REF] Cuypers | Virtual genomes in flux: an interplay of neutrality and adaptability explains genome expansion and streamlining[END_REF]]. Yet, EvoEvo is an integrative concept exactly as fitness is. Indeed fitness is the result of the interaction of all the organisation levels of the organism, including its interactions with its environment. Similarly, the robustness/evolvability/variability of the phenotype is the result of the interaction of robustness/evolvability/variability at all the organisation levels of the organism (including its interactions with its environment!). Furthermore, these properties are not independent and they may interact in a cooperative or competitive way (e.g. evolving chaperone proteins reduces the phenotypic variability, thus increasing the robustness). That is why a computational model of EvoEvo must be integrated, including the main organisation levels of the genotype-to-phenotype map (genome, transcription network, metabolic network, phenotype, fitness, population). Moreover, the genotype-tophenotype mapping must be evolvable, meaning that the complexity and struc-ture of all these levels must be able to change and that these changes must be able to induce changes of variability/robustness/evolvability. Such a model will necessarily incorporate a large set of parameters and its study is likely to be very difficult. However, it will give rise to new hypothesis and predictions, impossible to obtain with previous models. Indeed, such a model must be seen as a proof-of-concept model [START_REF] Servedio | Not Just a Theory?The Utility of Mathematical Models in Evolutionary Biology[END_REF], used to investigate EvoEvo theory and to test or to generate predictions.

Overview of the integrated evolutionary model

In the context of the EvoEvo project, we have designed an integrated model to study Evolution of Evolution. Obviously, such a model cannot include the whole complexity of real organisms. Moreover, one has to keep in mind that our objective is not to study the evolution of such or such organism. Its very aim is to study the evolutionary process and to unravel the EvoEvo strategies that result from a pure Darwininan evolution. That is why we propose to study evolution and EvoEvo in a simplified, abstract, world, designed by an "artificial chemistry" [START_REF] Dittrich | Artificial chemistries -a review[END_REF]. This artificial chemistry provides a set of objects and a set of rules that govern their interactions. In the model, organisms will be composed of these objects, the rules giving them their dynamic and, ultimately, their fitness.

To build our integrated model, we designed a modular artificial chemistry: the set of objects and the set of rules are split into modules that represent the different organisation levels we want to study as well as the interactions between those levels (e.g. a specific set of rules specifies how the genome is transcribed and translated into proteins). Here we propose to include five levels in our model: the genome, the genetic network, the metabolic network, the fitness and the environment (note that we don't include any "phenotype": the phenotype will simply be the result of the metabolic network dynamic in the organism's environment). These levels are described in the following sections.

Basically, our integrated evolutionary model is an individual based model. Each individual is an asexual virtual cell owning a genome. This genome encodes a regulatory network and a metabolic network. The metabolic network can uptake and convert nutrients from the environment, resulting in the ability (or not) for the cell to divide. Dividing cells form a population that grows on a two dimensional environmental grid providing fresh nutrients, but also nutrients and waste released by cells, actively or after death. This dynamical process results in a modification of the environment and allows for the emergence of complex ecosystems. It also creates the conditions for local competition between cells.

In the model, the molecular structure of the organisms is entirely defined by their genome (possibly with an interaction with the organism's environment). This genome can undergo mutation at each replication (point mutations and large rearrangements). The interaction between this variation process and the competition process described above results in a Darwinian evolution: individuals become more and more adapted to their environment and able to replicate more an more efficiently. However, by doing so, they also modify the shared environment (e.g. by releasing new metabolites), thus changing their own evolutionary conditions...

In this part, we present our formalism, inspired from [START_REF] Crombach | Evolution of evolvability in gene regulatory networks[END_REF] and [START_REF] Beslon | From digital genetics to knowledge discovery: Perspectives in genetic network understanding[END_REF], and describe the integrated evolutionary model. In the following section, we will present an example of the evolutionary dynamic observed in the model.

Description of the genome level

Genome structure. Each organism owns a circular string of functional/nonfunctional elements. The genome is a coarse-grained genome, inspired on [START_REF] Crombach | Evolution of evolvability in gene regulatory networks[END_REF], and defined as a list of functional or non-functional elements mathematically defined as n-tuples1 (the elements of our genome). Each functional tuple is parametrised to define its function and its number of dimensions n (depending on the type). Five types of tuples are defined in our artificial chemistry: E: Tuples coding for enzymes performing reactions in the metabolic network, via the following Michaelis-Menten equation:

d[p] dt = k cat .[E].[s] k m + [s] (1) with s, p ∈ N * , k cat ∈ R and k m ∈ R + being encoded in the tuple, [s] and [p]
the concentrations of the metabolites s and p, and [E] the enzymatic concentration (we assume that the concentration of free enzymes [E] is always equal to the total concentration [E T ]), TF: Tuples coding for transcription factors. Each transcription factor i owns a binding site identification tag j ∈ Z and an affinity A ij for this binding site, BS: Tuples coding for binding sites specifying which transcription factor may bind to them via their own identification tag ∈ Z, P: Tuples coding for promoters determining where the transcription should start. Each promoter i owns a basal expression level β i , NC: Non functional tuples constituting the non coding part of the genome.

Binding sites directly flanking a promoter regulate its transcriptional activity. The enhancer site directly precedes the promoter and is made of one or more contiguous binding sites. The operator site directly follows the promoter and is also made of one or more contiguous binding sites. Transcription factors that bind on the enhancer site increase the transcriptional activity. On the opposite, transcription factors that bind on the operator site down-regulate the promoter activity (see figure 1.1). As in R-aevol [START_REF] Beslon | From digital genetics to knowledge discovery: Perspectives in genetic network understanding[END_REF], a promoter has a basal level activity β. TF or E tuples following the operator site are transcribed, thereby allowing for operons. Downstream of the operator site, any tuple other than TF or E makes the transcription stop. To be functional, the promoter can be flanked by binding sites or not, but TF or E tuples must immediately follow the regulation unit (enhancer site + promoter + operator site, see figure 1.2). Mutational operators. The genome undergoes point mutations and large rearrangements during replication. If a tuple undergoes a point mutation, it operates a jump in the tuple space by adding a n-dimensional random vector. A tuple can be unfunctionalised by a point mutation or during a large rearrangement if it is located on a breakpoint. Non coding tuples can also be restored into one or another functional type, however it is impossible to mutate directly from a functional type to another (see figure 1.4). All the mutation rates are configurable. The genome also undergoes large chromosomic rearrangements: duplications, deletions, inversions, and translocations. The various types of mutation can modify existing tuples, but also create new tuples, delete some existing tuples, modify the length of the non coding regions, modify tuple order...

Description of the genetic regulation network level

The genetic regulation network (GRN) is computed from the interactions of transcription factors (TF) and binding sites (BS) elements. Its activity is computed in four steps:

1. The activity A s (t) of each binding site s is A s (t) = j c j (t).A js with c j (t) the concentration of the transcription factor j at time t and A js the affinity of this transcription factor for the binding site s, 2. We then compute the activity of the enhancer site E i (t) and of the operator site O i (t) flanking the promoter i:

   E i (t) = j∈enhanceri A s (t) O i (t) = j∈operatori A s (t) (2)
3. Then, the transcription rate e i over time of the promoter i is computed, possibly including a stochastic component [START_REF] Elowitz | Stochastic gene expression in a single cell[END_REF] represented by a transcriptional noise η genetically encoded in the promoter [START_REF] Roberts | Noise contributions in an inducible genetic switch: a whole-cell simulation study[END_REF]. Note that η mutates as all the elements of the tuple, allowing for evolution of stochastic gene expression [START_REF] Newman | Single-cell proteomic analysis of S. cerevisiae reveals the architecture of biological noise[END_REF]. Then the transcription rate e i of each promoter i is given by an Hill-like function:

e i (t) = β i . θ n O i (t) n + θ n . 1 + 1 β i -1 E i (t) n E i (t) n + θ n + ξ i (t) (3)
with β i the basal expression level of the promoter i, n and θ being constant coefficients that determine the shape of the Hill function. ξ i (t) is a random number drawn from the gaussian distribution N (0, η i ). 4. The transcription rate e i is applied to each E or TF tuple being controlled by the promoter i, such that each protein product (enzyme or transcription factor) has its own concentration regulated through a synthesis-degradation rule, depending on e i :

∂c i ∂t = e i (t) -φi(t) (4)
where φ is a temporal scaling constant representing the protein degradation rate.

Description of the metabolic network level

Enzymes performing reactions in the metabolic network are encoded by tuples of type E. If s = p, the reaction takes place in the cytoplasm of the cell. If s = p, the enzyme is a inflowing or outflowing pump depending on the sign of k cat . For each cell, the whole set of reactions defines an ordinary differential equations (ODE) system, which is solved numerically. Some metabolic products are essential for the cell's growth, and some other are intermediate products or waste. In the integrated evolutionary model, prime numbers are considered to be essential metabolites: their production contribute to the growth rate by increasing the probability to produce offspring. Overproducing metabolites can also lead to cell's toxicity. Hence, one can define toxicity thresholds for essential and non essential metabolites. Over-reaching a toxicity threshold impairs cell's fitness. Finally, during replication, daughter cells share cytoplasmic content at division (proteins and metabolites). It is also possible to define energy constraints in the artificial chemistry, such that cells must perform catabolic reactions to earn energy and produce essential metabolites.

Coupling the genetic and the metabolic networks

Bacteria are able to sense their environment by detecting the presence of a particular molecule or signal, and to give an appropriate answer by updating their gene expression profile. In the integrated evolutionary model, co-enzymes can repress or activate transcription factors activity. This is done by adding three elements to the transcription factor tuple (TF): A co-enzyme identification tag ∈ N * , a free activity (A f ree , boolean) and a bound activity (A bound , boolean).

A metabolite m ∈ N * can repress or activate a TF acting as a co-enzyme:

-If A bound = true and A f ree = f alse, the co-enzyme activates the TF, -If A f ree = true and A bound = f alse, the co-enzyme inhibits the TF, -If A f ree = f alse and A bound = f alse, the TF is never active, -If A f ree = true and A bound = true, the TF is always active.

Finally, the concentration c i (t) of the TF and the concentration coE i (t) of the co-enzyme are combined (depending on the values of A f ree and A bound ) to compute the active fraction of c i (t).

Description of the population and environment levels

Individuals "live" on a two dimensional grid, each grid site containing at most one individual. The physical environment is described at the grid level: each grid site contains a list of free metabolites, each with its concentration level. Those free metabolites diffuse with a diffusion rate D and are degraded with a degradation rate D g . Individuals compete for the free metabolites and to produce offspring in empty sites. Individuals interact with their local environment by pumping metabolites in and out and releasing their content at death. Metabolites can also diffuse through the cell membrane at rate D m . In this case, pumps are active mechanisms that the cell can use to maintain an internal concentration different from the external one.

At each simulation time step, organisms are evaluated and either killed, updated or replicated depending on their current state:

1. An active cell can go to death. Death probability follows a Poisson law. At death, cell content is released in the local environment, 2. If the cell do not die and is unable to divide (e.g. because there is no free space in its neighbourhood), its genetic regulation network and metabolic network are updated, and its score is computed from its metabolic concentrations state vector X through the score function F (e.g. the sum of essential metabolite concentrations). 3. For each empty grid site, all active neighbours compete to select the replicating individual, depending on relative fitnesses. To avoid biases, empty grid cells are updated in a random order.

This simple framework allows to model different real experimental setups, including serial plates or chemostat [START_REF] Mozhayskiy | Microbial evolution in vivo and in silico: methods and applications[END_REF]. Similarly, some individuals can be regularly picked up in the environment to seed a new colony, thus mimicking a mutation accumulation experiment.

Simulation example

The integrated evolutionary model allows the user to choose the complexity of the simulation by removing or including biological organisation levels. As an example, we present here the results of a simulation where cells evolved their metabolic network in a highly variable environment of size 32 * 32 for 100000 time steps (approximatively 10000 generations). Some features have been shut down: the genetic regulation network, the transcriptional noise, energy constraints and membrane permeability. The meaningful parameters of the simulation are displayed in table 1 (among 55 parameters). The score function F is the sum of essential metabolite concentrations produced by the cell, and is updated at each time step of the cell's life. Metabolite m = 1 at concentration [m] = 20.0 is introduced in the environment as an exogenous nutrient, at random sites with a probability p = 1e -05 per grid site per simulation time step, making the environment a low nutrient and variable one.

Results

The simulation provides a large set of outputs and statistics, available in a dedicated HTML viewer2 . Lineage and phylogenic trees are also dynamically computed during the simulation, and displayed in the viewer. Examining the Score. The score of each generation is averaged through the cell's lifespan, and is boxed by the minimum and maximum scores reached by the cell during its life. Here, the score tends to increase but is sometimes impaired (at approximatively 2000, 7000 and 10000 generations). Cytoplasm metabolic content. Each cell inherits an half of the cytoplasm of its mother (light blue), and then starts to accumulate nutrients (dark blue).

The score is strongly correlated with the evolution of the cytoplasm content. Genome size. Thanks to large rearrangements, the genome size varies through generations. Here, the genome undergoes strong size variations. Large genome size duplications seems to be correlated with score impairment. Proportions of functional (coding) and non-functional (non coding) elements indicate that non coding elements represent a large part of the genome. Genetic redundancy for metabolism. Several genes in the genome can code for the exact same enzyme in the metabolic network. This genetic redundancy associated with the metabolic network is strongly correlated with genome duplications and deletions.

To grow (i.e. increase their fitness), cells need to produce and accumulate essential metabolites. Yet, toxicity thresholds and limited nutrient concentrations in the environment strongly restrict fitness improvement. Thus, the best lineage progressively increases its fitness but sometimes undergoes rapid fitness drops (we call them "crisis"). To understand the origin of these events, we now examine population and environment levels. The population size is variable, alternating plateauing phases not far from the maximum grid capacity (e.g. between 35000 and 50000 time steps), and size falls during crises (e.g. at 60000 or 100000 time steps, see figure 3.1). Evolution of environmental metabolic concentrations shows a rapid diversification (see figure 3.2), with temporary accumulation of one specific essential metabolite at high levels, followed by metabolic sweeps. Sweeping periods are associated with population crisis. To investigate more deeply simulation results, the integrated evolutionary model allows us to go back in the evolutionary history by "reviving" ancestor populations saved in backup files. Each 1000 time steps, we computed the trophic network of the population, classifying cells depending on their interaction with other cells and their environment. The cell's state (cytoplasm content, enzymatic reactions, inflowing and outflowing pumps) defines a "profile" used to build a graph, each node corresponding to one specific profile (or group). Edges of the graph represent consuming (solid edges) or helping (dashed edges) relationships. We distinguish three trophic levels: Level 0: Trophic groups eating and converting nutrients that are not produced by other individuals. These nutrients can be m = 1 (the nutrient provided randomly) or any other nutrient released previously by a now extinct colony, Level 1: Trophic groups eating environment nutrients and the product of other cells, Level 2: Trophic groups only eating nutrients actively produced by other cells.

We recovered the proportion of groups belonging to each level (figure 3.3) and the numbers of cells per group (figure 3.4). We also recovered the trophic level of the last best individual lineage, as shown at the bottom of the figure 3 indicates that the number of groups per level fluctuates, increasing during crisis, and decreasing in between. The number of groups belonging to level 0 tends to be maintained at low levels. Between crisis, the level 2 is the most representative one in the population, but the situation is reverted during crisis, where levels 0 and 1 increase (see figure 3 [START_REF] Wilke | Evolution of digital organisms at high mutation rates leads to survival of the flattest[END_REF]. The trophic network 3.5.3 (at 95000 time steps, during a plateauing phase) shows a complexity reduction, but still maintained at a higher level than for the trophic network 3.5.1.

Discussion

At the beginning of the simulation, cells quickly evolve a metabolic network that uptakes the exogenous nutrient (m = 1, in light brown on figure 3.2), and convert it in essential metabolites. The task is tough, the nutrient being rare and localised, so the population is at low level. At death, cells release their cytoplasm content in the environment, thus enriching it. This accumulation of nutrients is followed by the rapid growth of level 2 cells, specialised in consuming the product of others cells. The decline of level 0 and level 1 cells leads to the depletion of level 2 cells nutrient, and the whole population decreases. At some point, new level 0 and level 1 cells convert the exogenous metabolite in a new essential one, cycling this behavior. Some questions arise from those results:

(1) Why don't we observe the fixation of level 1 cells by niche exclusion? The reason could be the cost to convert the exogenous metabolite, such that level 2 cells always grow faster if nutrients are available,

(2) Why do we observe a innovation process at each population crisis, instead of a simple negative frequency-dependent behaviour? If level 0 and level 1 cells compete with level 2 ones for available space, a strategy could be to innovate and produce unseen metabolites that level 2 cells cannot consume. The evolution of concentrations in the environment (figure 3.2) could support this conjecture. Indeed, new essential metabolite production always starts at the very beginning of population crisis,

(3) What is the role of large genome size duplications during crisis ? And what is the link with trophic network complexity variations ? Duplicating the genome size increases enzymatic dosage (enzymatic reactions go faster) but also the mutation rates. If a competition occurs between levels, indirect selection could favour mutators that try to escape their niche to find a new metabolic strategy.

These questions are clearly relevant in an EvoEvo perspective and this preliminary study can be continued by e.g. comparing these results with simulation in which the same amount of food is provided but in an homogeneous way. Similarly, the effect of genomes duplication must be investigated further to understand more precisely their origin and consequences.

Perspectives

In this example, we observed the emergence of a complex negative frequencydependent behaviour, associated with regular metabolic innovations. To understand it, we investigated several organisation levels of the model, linking for example the evolution of the genome structure to the metabolic network and the ecosystem levels. Our preliminary results remind situations observed in bacteria, e.g. the negative frequency-dependent strategy observed in one of the Lenski E. coli lineages [START_REF] Rozen | Death and cannibalism in a seasonal environment facilitate bacterial coexistence[END_REF].

The integrated evolutionary model offers complex evolutionary outcomes, but also the tools to decipher them (regular backups, phylogenies, a large set of statistics, the recovering of trophic networks...). The complexity of the genotypeto-phenotype mapping and of the fitness landscape allows us to investigate the evolutive interactions of a large set of biological organisation levels, and then to deeply explore EvoEvo. Several questions could be investigated, such that the evolution of robustness, evolvability or open-endedness, with the hope to open the door towards a unified theory of Evolution of Evolution.

Fig. 1 . ( 1 )

 11 Fig. 1. (1)A functional region starts with a promoter, possibly flanked by an enhancer site and/or an operator site. All contiguous E or TF tuples following the operator site are transcribed. The first tuple of another type interrupts the transcription (here a NC tuple). In this case, the functional region is an operon. (2.A) The promoter can be flanked by binding sites or not, but TF or E tuples must immediately follow the regulation unit (enhancer + promoter + operator). (2.B) On left, a non coding tuple interrupts the transcription. On right, the promoter is missing. (3) The genome is a circular single-strand sequence of tuples. At each replication, the genome undergoes mutations: (A) large duplications, (B) large deletions, (C) translocations, (D) inversions. Red arrows symbolize breakpoints in the sequence. (E) point mutations and breakpoints can unfunctionalise or functionalise a tuple. (4) Functional types (E, TF, BS and P) can be unfunctionalised with probability p NC . A NC tuple can be restored to one type or another depending on 4 mutations rates: p BS is the probability to be restored in a binding site (resp. p PROM , p E and p TF ).

Fig. 2 .

 2 Fig. 2. Some variables showing the evolution of the last best individual lineage. 1. Evolution of the cell's score. 2. Evolution of the cytoplasm metabolic content. 3. Evolution of the genome size. 4. Evolution of the genetic redundancy associated to the metabolism. The meaning of each variable is explained in the results section 3.1.

  .

4 .

 4 The three levels undergo strong fluctuations during the simulation.

Figure 3 Fig. 3 .

 33 Fig. 3. Evolution of the population and the environment. 1. Evolution of the population size. 2. Evolution of metabolic concentrations in the environment. 3. Evolution of the number of trophic groups per trophic level. 4. Evolution of the number of cells per trophic level. 5. Representation of three typical trophic networks. Node labels represent metabolites pumped in by the trophic group (the "uptake profile"). Solid edges represent a consuming activity (e.g. in 5.1, group 19 consumes group 7 metabolic products). Dashed edges represent "helping" activity, i.e. active release. "-" groups have no inflowing pumps. "ENV" group is the environment, taking only into account the exogenous nutrient.

Table 1 .

 1 Meaningful parameters used for the simulation detailed in this section, among 55 parameters. lineage of the last best individual of the simulation gives a first idea of the evolution of main simulation variables, as exemplified in the figure2:

	Parameter	Symbol Value Parameter	Symbol Value
	Simulation time		1e+05 Stochastic gene expression		NO
	Energy constraints		NO Membrane permeability		NO
	Metabolic inheritance		YES Selection threshold	θ	0.6
	Selection pressure	p	10 Grid width	W	32
	Point mutation rate		1e-04 Grid height	H	32
	Duplication rate		1e-04 Unfunctionalisation rate	pNC 1e-04
	Deletion rate		1e-04 Functionalisation rate	pE	1e-06
	Translocation rate		1e-04 Metabolites toxicity threshold		1
	Inversion rate		1e-04 Death probability	p death 0.05
	Diffusion rate	D	0.01 Degradation rate	Dg	0.001

  .4). The bottom part of figure3.4 indicates that trophic groups are not monophyletic. Indeed, the lineage of the last best individual regularly transits between levels, even if level 2 is the most frequent one (grey level represents cells consuming environment wastes, not produced anymore by other cells). Variations in figure3.3 is supported by trophic network examples. In figure3.5.1 (at 40000 time steps, during a plateauing phase), interacting groups are only level 2 groups. Level 0 groups do not produce anymore nutrients for those groups (m = 7 and m = 19), that are thus devoted to extinction. The trophic network 3.5.2 (at 62000 time steps, during a crisis) shows a large complexification of the network, with a new metabolite being produced (m =

A n-tuple is an ordered list (x1, x2, ..., xn) : T1 × T2 × ... × Tn with Ti the product type of xi (e.g. R, N, ...).

The viewer of this simulation is available at http://liris.cnrs.fr/ ~crocaber/ SimulationViewer/viewer/viewer.html

Acknowledgments This work was supported by the European Commission 7th Framework Programme (FP7-ICT-2013.9.6 FET Proactive: Evolving Living Technologies) EvoEvo project (ICT-610427).


Gérard

« Effets de ségrégation dans les procédés utilisant des mélanges de poudres » « Segregation effects in processes handling powder mixtures »

----------------------------------------------------------------------------------------------------------------

Conferences Chem and Ind

Paris, 13-15 septembre 2005, Cité des Sciences et de l'Industrie (webs_limousin@yahoo.fr)

----------------------------------------------------------------------------------------------------------------

The handling of particles is a general operation concerning a great number of industries for which a process implies a treatment of matter or material: foods, pharmaceutical, cosmetics, chemical. This represents 60-70% of the processes. The materials can be of any type, mineral, organic.

The handling of materials is just a transfer from a place to another one, thanks to different devices by operating under different conditions.

This insures the transfer of particles from trucks to a bin or from a bin to the unit where the transformation of particles will occur, for example from mixer to bin, or bin to die where compaction or tabletting will occur, or bin to bags.

Generally speaking a process can be described through a series of operation units allowing the suitable transformations, like elaboration of particles by crystallisation, heating, melting, sintering or formulation by mixing of particles, leading to the material presenting the properties satisfactory for the end-user.

Unfortunately handling particles is not considered as a unit operation: nothing is intentionally done to change the material so that one supposes that nothing bad can happen. In particular the storage of particles that could appear as a totally neutral operation may be completely neglected by engineers searching causes of a variability in the properties of a product.

Nevertheless the powders handled in industries present quite different behaviours, depending on the characteristics of the particles themselves or the conditions under which they are used. In particular the segregation of particles, i.e. the spatial separation of 18/01/07 1 GT/ENSM-SE/SPIN/LPMG particles making the granular medium composition not homogeneous is a phenomenon largely encountered.

A review of the main mechanisms of segregation will be given and the effects of both parameters -particle characteristics and conditions-upon the segregation of particles will be examined. Some examples from laboratories experiments and industrial problems illustrate the difficulties encountered by engineers to control and avoid serious inconveniences occurring in processes.

Among the different parameters characterising the particles, the one exhibiting the strongest influence in this kind of problem is the size, represented by the radius R of the sphere equivalent to the particle. The particle size distribution is one of the keys for the knowledge. A second is the adhesion energy Wa between a wall and a particle or the cohesion energy Wc between two particles. The particle shape and density can be also related to segregation effects. Of such effects, size difference is by far the most important, and almost always this is the most serious cause of segregation. In comparison the other effects are minor and unimportant [START_REF]Mixing in the process industries[END_REF].

The key parameters for the external parameters defining the storage conditions are mainly the temperature or relative humidity controlling the surface properties, and adhesion or cohesion. The geometrical parameters defining the shape and dimensions of the storage devices like silos are very often calculated in order to optimize the solution to flow problems: arching, piping, or segregation due to particle size. In addition mechanical energy can be also considered to facilitate inlet or outlet operations. For these reasons, the studies must take into account two main types of powders: cohesive or free flowing particles [START_REF]Handbook of powder technology[END_REF]. A given powder can exhibit a cohesive behaviour whenever its radius is lower than a critical one Rc, or a free flowing behaviour for R>Rc. The critical radius Rc depends on the particle interactions and the particle density.

A great number of experiments and discussions have been devoted to segregation effects in particle mixtures. Anyway the understanding of segregation is not achieved at all, even if a series of models allows some guide for interpretation.

Hereafter let us consider the main mechanisms of segregation.

A first type of segregation appears in particulate systems where the particles can behave without grain-grain interactions, i.e. in lean phase systems. The segregation by projection, elutriation, or decantation is the result of the forces exerted on the particlegravity, drag, friction and driving forces.

This may yield opposite results: in air without a driving stream imposed to a particle, a small particle will be stopped very quickly after being projected from a gun. Whereas when a laminar gas flow acts continuously on small and large particles, the small particle will go further than a large one. Such size dependent movements can be also at the basis of particle size determination: sedimentation, fluorometry in lab, gravitational separation units in industry.

Other segregations take place in dense phase systems for packing with a packing fraction greater than 0.2.

the segregation by sifting corresponds to screen effects and may be observed when pouring a heap. Large particles run down to the edge of the heap and small particles remain in the middle: small particles are more easily stopped when going down the heap because big particles offer a higher barrier and trapping sites for small particles, and the rotational or translation kinetic energy of small particles is not enough to let them go further.

In the shear segregation, the percolation of small particles through a granular packing is created by shear stresses imposed to the bed. In this case, some mobility can be bring to the particles by an external device (drum shaker for example) and if a velocity gradient is set up any one layer of particles have a velocity relative to a neighbouring layer. If an opportunity is presented, particles can roll or percolate into a new layer. In fact smaller particles have a greater mobility because they need smaller sites in the adjacent layer, and segregation can take place. Such velocity gradients can be seen in milling or mixing operations in drums: 2D experiments show that small particles appear in the centre band of the bed and large particles stay around them (radial segregation). 3D experiments can also give sometimes another spatial separation, in the axial direction of the cylindrical drum, with alternate bands of the different kinds of particles.

In the segregation due to vibration, percolation of small particles downside through a granular packing can undergo also almost spontaneously if the ratio of the large on small particles is greater than 6, or it may be created by vibrations. This can be evidenced by the experiment in which a large and dense ball is placed at the bottom of smaller sand particles: after vibrations or shocks it will rise to the top surface after a few seconds.

The segregation effects can be characterised by two indexes: the intensity and the scale of segregation defined as well in mixing studies. The first expresses the composition variation from a segregated region to a region having the mean composition; the second indicates the length -or volume-separating the domains unmixed, each having a given composition.

They mainly depend on the powder properties that can be precised through laboratory tests like volumenometry, angle of repose, shear cell, image analysis, granulometry… Of course these effects of segregation and transformation of particle-particle interactions can modify the composition of the granular media before, during or after the storage and transfer operations.

For storage or feeding applications two kinds of systems can be considered: lean phase flow or dense particle flows, corresponding respectively to industrial pneumatic conveying of particles dispersed in a gas flow of high velocity -with small concentration of solids-, or screw conveyors, belt feeders, volumetric proportioning equipments, plug conveyance -dense phase conveying systems-. As seen before, some relevant difficulties can be encountered linked to segregation effects.

In silos or hoppers for example, the user is always in front of a squeeze: to discharge a silo, discharge aids are intended to improve the bulk solid flow and draw the bulk solid out of the discharge outlet. For that purpose, pneumatic or mechanical devices (air injection, stirrers), vibrations by shakers or flow agents can be added. In the same time possible segregations might appear.

Vibrating conveyors are not suited for easily fluidisable solids.

The storage in a bin can be affected by the vicinity of external sources of vibrations.

18/01/07