Aimé Lachapelle

Marie-Therese Wolfram

On a mean field game approach modeling congestion and aversion in pedestrian crowds

Keywords: Mean field games, interacting populations, Nash equilibrium, rational expectations, flow of pedestrians, lane formation, numerical approximation

In this paper we present a new class of pedestrian crowd models based on the mean field games theory introduced by Lasry and Lions in 2006. This macroscopic approach is based on a microscopic model, that considers smart pedestrians who rationally interact and anticipate the future. This leads to a forward-backward structure in time. We focus on two-population interactions and validate the modeling with simple examples such as self-organization behavior as for instance lane formation. Two complementary classes of problems are addressed, namely the case of crowd aversion and the one of congestion. In both cases we describe the model, build a numerical solver (respectively based on optimization formulation and partial differential equations), and finally provide some numerical tests involving complex group behaviors such as symmetry breaking and lane formation.

Introduction

In the last decades crowd motion modeling has become an active area of research. As the world population and the urbanization continues to increase we can not underestimate the importance of understanding the behavior of human crowds. Nowadays there is lot of interdisciplinary research among physicists, sociologists, biologists and mathematicians; crowd motion has become one of the emerging research topics. The first empirical studies on crowd motion have been deducted about 50 years ago (cf. [START_REF] Hankin | Passenger flow in subways[END_REF]; [START_REF] Hoel | Pedestrian travel rates in central business districts[END_REF]), nowadays various modeling approaches can be found in the literature. We distinguish between two different approaches, either microscopic models or macroscopic ones. The most popular microscopic approaches include the behavioral force models (cf. [START_REF] Helbing | Microscopic foundation of stochastic game dynamical equations[END_REF][START_REF] Helbing | A mathematical model for the behavior of pedestrians[END_REF]; [START_REF] Helbing | Simulation of pedestrian crowds in normal and evacuation situations[END_REF]; [START_REF] Helbing | On the fluid mechanics of human crowd motion[END_REF]), cellular automata (cf. [START_REF] Burstedde | Cellular automaton approach to pedestrian dynamics-applications[END_REF]; [START_REF] Kirchner | Cellular automaton simulations of pedestrian dynamics and evacuation processes[END_REF]) and stochastic models. Behavioral force models have been very successful in describing various crowd motion phenomena like lane formation, oscillations at bottlenecks or clogging. The main drawback of microscopic approaches is the high computational cost for large number of particles. In macroscopic models the crowd is treated as a density. Here the computational cost is much smaller, but it's less clear how to describe the motion of the density correctly. Several macroscopic models have been derived by considering the limit as the number of particles tends to infinity, e.g. Burger et al. (submitted). Other macroscopic approaches are based on kinetic equations or fluid dynamics (see e.g. [START_REF] Helbing | On the fluid mechanics of human crowd motion[END_REF]). A different macroscopic approach has been proposed by [START_REF] Hughes | A continuum theory for the flow of pedestrians[END_REF], [START_REF] Francesco | On the hughes' model of pedestrian flow: The one-dimensional case[END_REF] studied the analysis and numerics in the one-dimensional case. More recent models are based on optimal transportation, see Maury et al. (to appear). For a general overview on the various modeling approaches we refer to [START_REF] Helbing | Simulation of pedestrian crowds in normal and evacuation situations[END_REF]. Furthermore we would like to point out the connection between human crowd motion models and herding models for animals, like ants, sheep or fish. A model for describing the collective motion of ants has been used successfully to describe lane formation in crowds. Ants use chemical substrates to lead others to valuable food resources. [START_REF] Schadschneider | Pedestrian and Evacuation Dynamics[END_REF] used a cellular automata model which uses a virtual trace that features collective behavior of humans an subsequent lane formation.

In this paper we would like to present a new approach to model crowd motions. It is called mean field games (MFG) and has been introduced initially by Lasry and Lions, see [START_REF] Lasry | Mean field games[END_REF]. A particular appealing feature is that the approximating macroscopic mean field model is derived from a mean-field game as the number of players tends to infinity. Consequently it offers the tractability of a macroscopic framework together with the more realistic interpretation at the microscopic level. Furthermore we treat players as "smart" individuals which try to optimize their path with respect to a particular goal. Thus, a MFG model is quite novel for crowd modeling and differentiates itself from the models cited above in mainly two points. First we don't talk anymore about neither particles nor robots nor automata, but rather about real individuals or agents, having strategical interactions within the crowd. Second, these individuals anticipate the future. This is mathematically expressed through a forward-backward structure, the forward dynamic describing the crowd dynamic whereas the backward one is needed to build the expectations. This differs from pure forward frameworks that might be unable to describe particular features of crowd dynamics. We believe that a MFG model is interesting in the sense that beyond offering a description of how pedestrian behave in a crowd, it also answers the questions why and how the crowd moves (or why it has its specific shape in the stationary case). Here is an example of why it could be useful: if for any exogenous reason some people in the crowd suddenly stop for a while, the collective behavior will be re-established. This cannot happen in automata models. As indicated in the name mean field games, it is a game approach and has the very important feature to consider the pedestrians' behavior as a non-cooperative equilibrium. Therefore it is unfortunately not possible to consider cases of small groups interacting within the crowd.

We reiterate that the limiting MFG are motivated by an N-player stochastic game. These limiting macroscopic partial differential equations (PDEs) have a surprisingly simple structure, which allow the development of efficient numerical methods. Different numerical approaches for MFG can be found in the literature (see for instance [START_REF] Achdou | Mean field games: Numerical methods[END_REF]; [START_REF] Guéant | A reference case for mean field games[END_REF]; [START_REF] Lachapelle | Control and transport problems in economics: Theoretical and numerical aspects[END_REF]; [START_REF] Lachapelle | Computation of mean field equilibria in economics[END_REF]; Lions (2008Lions ( -2010))). They are either based on the optimal control formulation of the problem or on the PDE system. In the first case different methods from parabolic optimal control theory are available. On the PDE level Newton-type methods have been used successfully for numerical simulations. Nevertheless the forward-backward structure and the highly nonlinear equations still pose a great challenge to numerical analysts and the development of efficient numerical methods is still in its early stages.

This paper is organized as follows. Section 2 is devoted to the general introduction of mean field games in a crowd dynamic framework. Next we turn to the central topic of our paper, that is the study of two interacting populations (or equivalently crowds). More particularly we present two different classes of MFG (in order to show the diversity of cases that could be treated by MFG). For both of them we describe the modeling, discuss a numerical solver and provide numerical results. In section 3 we present a model where individuals are averse to the crowd, in a time-dependent framework. The gradient-descent numerical method is based on an optimal control formulation of the problem. Our numerical tests mainly provide a symmetry breaking in some xenophobic situation (where the pedestrians of the first group really dislike the pedestrians of the other group). In section 4 we look at a two-species stationary model involving congestion.

Here the numerical simulations, based on Newtons' method, reproduce well known phenomena in crowd motion like lane formation. Finally, we close with some concluding remarks in section 5.

Mean field games: a micro-macro approach to model pedestrian crowds

In this section we discuss the general concept of mean field games and its application to crowd motion dynamics. For reasons of clarity and readability we start with the introduction of mean field games for single species first. This discussion can easily be generalized for two populations, which is of particular interest when studying interactions among crowds. First we recall the microscopic fundament of MFG, that is a multi-player stochastic differential games. Then we present the approximation of such games when the number of player tends to infinity (without entering the mathematical details), namely the associated mean field game. Various remarks throughout the section have been added to clarify the modeling assumptions.

Microscopic fundament: stochastic differential games

Let us start at the microscopic scale of a MFG, namely an N-player stochastic differential game. We look at an evolution game on a period [0, T ] in a stochastic framework. Here we choose the probability space (Ω, F, P) where F = (F t ) t∈[0,T ] is the filtration enhanced by a vector (with length

N) of d-dimensional independent Brownian motions W = (W t ) t∈[0,T ] = (W 1 t ) t∈[0,T ] , ...(W N t ) t∈[0,T ] .
We consider N players (also referred as pedestrians or individuals) interacting one with each other through the choice of an action (or control). At every time t player i chooses an action α i from the set of feasible actions A i ⊂ R d , for i = 1, ..., N. We define A = A 1 × ... × A N and the global set of admissible actions as

A = {α = (α t ) t∈[0,T ] ; α is progressively measurable in A}.
Consequently, (α i t ) t∈[0,T ] denotes the strategy or control process of player i. For any α = (α 1 , ..., α N ) we will use the very convenient and usual notation in game theory (α i , α -i ) = (α 1 , ..., α i-1 , α i , α i+1 , ..., α N ). We are interested in simple cases where the state evolution of the system in R N×d is given by the solution of the following Stochastic Differential Equation (SDE for short):

dX t = α t dt + σdW t , X 0 = x ∈ R N×d , ( 1
)
where σ is a Nd × Nd diagonal matrix (independence case).

Remark 1. In other words player i evolves controlling the drift of the controlled stochastic process dX i t = α i t dt + σdW i t , starting at x i . Note that we will focus on the realistic particular case where the crowd evolves in a bounded subset of R 2 (in this case we will specify the boundary conditions later). We reiterate that the speed of each individual is not fixed a priori but rather determined by each pedestrian as the result of an optimization process (it is the control α). We come back on this comparison later, when we will describe the macroscopic setting. We emphasize that pedestrians are subject to independent noise. If we would consider correlated noises, the resulting mean field game would be much more complicated, see e.g. Lions (2008Lions ( -2010)).

It is well-known that for every α ∈ A there exists a unique solution X t = X α t to (1). Next we discuss how player i determines its optimal strategy by a cost functional. This cost functional consists of two parts, the first one corresponds to the running cost

f i : [0, T ] × Ω × R N×d × A i → R
and the second one to the terminal cost g i : Ω × R N×d → R. Then for each α ∈ A we define the cost functional

J i (α) = E T 0 f i (t, X α t , α i t )dt + g i (X α T ) . (2)
Note that we will also study the stationary case for which we will consider a slightly different definition of the cost functional namely:

J s i (α) = lim T →∞ 1 T E T 0 f i (t, X α t , α i t )dt .
(3)

Next we discuss the special case of mean-field type interactions between players i.e. we discuss special forms of costs and their very important aspect in the modeling.

Each pedestrian takes into account her own state X i t and the mean-field created by the others, i.e. the empirical distribution of pedestrians x-i = 1 N-1 j i δ x j where δ x denotes the Dirac mass located at x. Then for each i and x = (x 1 , ..., x N ), the running and terminal costs have the particular form f i (t, x, α) = f (t, x i , x-i , α) and g i (x) = g(x i , x-i ) and we have

J i (α) = E T 0 f (t, X i,α t , X-i,α t , α i t )dt + g(X i,α T , X-i,α t ) .
Remark 2. Here we make the implicit assumptions that each player is anonymous. Indeed they are identical but can be different (that is we cannot identify who's who in the crowd but everyone can have her own path that may differ from the others). This can be seen through the form of the costs. For a fixed index i, J i is invariant by any permutation of the -i other players.

In our model we assume that the pedestrians are rational, with rational expectations, i.e. every individual or agent knows and anticipates the state of everybody through the distribution of the crowd. These assumptions are probably wrong in the real world (except in routine crowd situations, say the daily way towards the subway exit). Indeed the pedestrians could be partially and insufficiently informed, or some irrational local behavior can happen (e.g. an individual suddenly stops to observe the crowd). But as it is the case in physics, we first have to introduce the model with rational expectations and rational individuals in order to be able to correct and improve the model afterwards (e.g.

by introducing partial blindness or an extra noise for instance). A robot model might be more intuitive at the first glance, but it might fail if the environment changes.

A first step to consider partial blindness can be performed either through the introduction of a discount factor (chosen significantly high) which has been done in section 4 or considering running costs involving a convolution. The convolution of the empirical density with the indicator function of a circle of radius r corresponds to taking actions of individuals in certain neighborhood into account. Indeed we do not enter the details of the dependency of the cost upon the crowd empirical distribution, we do not even specified whether it is local or non-local. Here we will only consider fully non-local dependencies, the generalization to other cost functions will be subject of further research.

We reiterate that mean field games are an equilibrium-based approach. Hence we try to find a Nash equilibrium (or Nash point). A Nash point is a situation in which every player cannot improve his payoff by any unilateral deviation. That is a strategy α such that

∀i ∈ {1, ...N}, ∀α i ∈ A i , J i (α i , α -i ) ≥ J i (α).
Existence of Nash equilibria can be shown under very general assumptions on the cost. We will only consider quadratic dependence of f i upon α i , which is covered by standard existence results. A rigorous study of such problems can be found in [START_REF] Bensoussan | Nonlinear elliptic systems in stochastic game theory[END_REF]. Furthermore they show that the differential game is equivalent to a system of coupled PDEs, that are intractable when N is large. Next we discuss an approximation of such games by a mean field game with a continuous distribution of players.

Macroscopic framework: Mean Field Game equations

Though the notion of games with a continuum of players is not novel (see the work of Aumann ( 1964)), it has not been adapted to stochastic games (the difficulty comes from the independence of white noises when taking N → ∞). Our aim is not to enter the mathematical details of taking limits when N → ∞, we rather give a common intuition of how to derive the MFG system in the continuum of player setting. For more details on the highly technical proof of the approximation result, we refer to [START_REF] Lasry | Mean field games[END_REF] and Lions (2008Lions ( -2010)). Therefore we directly jump to the continuum of players setting. First we present the equations of the MFG system (whose solutions are mean-field equilibria or equivalently approximation of Nash equilibria) both in the dynamical and stationary cases. Then we briefly motivate a common intuition on how to derive MFG equations.

Dynamical and stationary MFG systems

We distinguish between finite horizon MFG and stationary MFG. This is motivated by the two differential crowd motion problems which we study in this paper. In section 3 we look at crowd motion dynamics and therefore consider the following (and more general) type of finite horizon MFG system:

∂ t v + σ 2 2 ∆v + H(t, x, ∇v, m) = 0, v| t=T = g(m T ) (4) ∂ t m - σ 2 2 ∆m + div(m ∂ p H(t, x, ∇v, m)) = 0, m| t=0 = m 0 . (5)
The mathematical structure of this system captures many features of MFG modeling. It is a forward-backward system coupling two PDEs. Equation ( 4) is a Hamilton-Jacobi-Bellman (HJB) equation on the adjoint variable v. It is related to the control aspects and basically gives the optimal strategy obtained by a backward reasoning (feedback). Note that there are two sources of coupling in this equation: both in the Hamiltonian H and in the terminal condition of the adjoint variable v. Here the Hamiltonian H is defined as the Legendre transform of the running cost f , i.e.

H(t, x, p, m) = f * (t, x, p, m) := sup α {pα -f (t, x, α, m)}.
On the other hand, equation ( 5) describes the mass evolution equation of the system. It is called a Kolmogorov equation in financial mathematics or a Fokker-Planck equation in physics. Here m 0 denotes the initial pedestrian distribution. Its solution m t is the distribution of pedestrians transported according to individuals' optimal choices (forward in time). The optimal strategy is given by α = ∂ p H(., ∇v, m).

The forward-backward structure is a core aspect of MFG. The achievement of a MFG equilibrium (i.e. a solution of the system (4)-( 5)) is performed through the rational expectations assumption (we come back more in details on this very important point in the next subsection).

In section 4 we look at the corresponding stationary case that is a MFG system having the following form:

σ 2 2 ∆v + H(x, ∇v, m) -rv + λ = 0 (6) - σ 2 2 ∆m + div(m ∂ p H(x, ∇v, m)) = 0 (7) λ ∈ R, mdx = 1, m > 0, vdx = 0. ( 8
)
The structure of the stationary problem is very similar to the the time-dependent one. Note the additional variable λ, which gives the optimal value of the equivalent of (3) in the game with a continuum of players, and the distribution's conditions in (8). The linear term rv appearing in (6) comes from the introduction of a discount factor r in the cost functional. We reiterate that the rigorous derivation of these PDEs systems can be found in [START_REF] Lasry | Mean field games[END_REF] and Lions (2008Lions ( -2010)).

Formal derivation of the finite horizon system

In order to better understand the forward-backward structure, we would like to give an intuition on the derivation of MFG. We consider an equivalent formulation of the individual problem consisting of minimizing (2) in the infinite number of players game. Here the starting point is the individual problem of a single individual located at x:

       inf α∈A E T 0 f (t, X α,x t , α t , m t )dt + g[m T ](X α,x T ) dX α,x t = α t dt + σdW t , X x 0 = x. ( 9
)
Note that we replace the empirical distribution 1

N-1 j i δ X j t by m t . Each individual in the continuum takes into account the global distribution evolution (m t ) t∈[0,T ] , hence anticipating the future. We proceed in two steps in order to get the system (4)-( 5). First we assume that (m t ) t∈[0,T ] is known (expectations hypothesis) and solve the stochastic control problem. We introduce the value function:

v(t, x) = inf α∈A E T t f (s, X α,x s , α s , m s )ds + g[m T ](X α,x T ) .
It is well known that the value function v is a solution to the HJB equation (see [START_REF] Bensoussan | Nonlinear elliptic systems in stochastic game theory[END_REF])

∂ t v + σ 2 2 ∆v + H(t, x, ∇v, m) = 0 , v| t=T = g[m T ].
We reiterate that we define the optimal feedback as a maximizer in the previous formula, computed at the point (t, x, ∇v(t, x), m t (x)), i.e. α(t, x, m) = ∂ p H(t, x, ∇v(t, x), m t (x)). The Kolmogorov equation describes the evolution of m t , it is a transport by drift (here α) and diffusion equation. For each t ∈ [0, T ], define the distribution as follows:

ϕdm t = E ϕ(X α,x t ) dm 0 (x) for all test function ϕ.
Then the Kolmogorov equation follows from applying Itô's Lemma

dϕ(X x t ) = [α(t, X x t , m t )∇ϕ(X x t ) + σ 2 2 ∆ϕ(X x t )]dt + σ∇ϕ(X x t )dW t .
Indeed, by a simple integration by parts one can write

T 0 ϕ∂ t m = d dt T 0 ϕdm t = T 0 E[dϕ(X x t )]dm 0 (x) = T 0 E[α(t, X x t )∇ϕ(X x t ) + σ 2 2 ∆ϕ(X x t )]dm 0 (x) = T 0 ϕ[-div(mα) + σ 2 2 ∆m],
and we finally obtain the desired equation. Now we are able to put everything together. For α(t, x) = ∂ p H(x, ∇v(t, x), m t (x)) we obtain the Kolmogorov equation:

∂ t m - σ 2 2 ∆m + div(m ∂ p H(t, x, ∇v, m)) = 0, m| t=0 = m 0 .
Remark 3. The forward-backward aspect can be considered as the heart of the MFG theory and refers to the rational expectations equilibrium approach. Indeed individuals anticipate the crowd evolution first and then evaluate their cost function ( 9). Next they deduce their strategy (feedback control) by solving the HJB equation. Finally the mass evolves according to these strategies (Kolmogorov equation). At the optimum the mass evolution has to coincide with the one which has been anticipated, according to the rational expectations assumption. Even in very simple configurations, a pure forward model should not be able to predict the crowd dynamic. This is why we believe that it is necessary to introduce forward-backward models such as MFG. Indeed, if we want to find expectation phenomenon, it is necessary to have both a forward equation (to describe the real dynamic) and a backward equation (to build up the expectations).

Link with optimal control

Finally we would like to complete the picture by linking mean field games to optimal control problems. MFGs are equivalent to optimal control problem under certain assumptions on the cost function. One of the presented crowd motion models is equivalent to an optimal control problem, the other not. In the first model, namely the crowd aversion model of section 3, we take full advantage of the optimal control structure in the numerical discretization. However, this won't be possible in the second model of section 4, where we take into account congestion.

An optimal control problem is equivalent to a MFG if the following criteria are satisfied. First the running cost can be written f (t, x, α, m) = L(t, x, α)m + V[x, m] and second V and g are the derivatives (e.g. the Gâteaux derivatives) of two potentials Φ and Ψ on bounded measures, i.e. V = Φ ′ and g = Ψ ′ . In this case, the critical points of the optimal control problem of Kolmogorov equation

       inf α∈A J(α) := T 0 Ω L(x, α)m(t, x)dx + Φ(m)(t)dt + Ψ(m T ) ∂ t m -σ 2 2 ∆m + div(αm) = 0 , m(0, .) = m 0 (.), (10)
are solutions of the MFG system (4)-( 5). A proof of this statement can be found in [START_REF] Lachapelle | Control and transport problems in economics: Theoretical and numerical aspects[END_REF]. Note that this is a sufficient condition as soon as Φ and Ψ are convex. Consequently, we see that the continuum of individual problems can be seen as a global optimization (or control) problem.

Finally we present some concluding remarks on our modeling approach.

Remark 4. First, it is clear that when we use the terminology of game we refer to classical game theory with an equilibrium approach and not to evolutionary game theory used in biology (also for crowds modeling and socio-behavioral modeling, see [START_REF] Haag | Concepts and models of a quantitative sociology[END_REF] or [START_REF] Helbing | Microscopic foundation of stochastic game dynamical equations[END_REF]).

• We reiterate one more time the difference with social force models (see [START_REF] Helbing | A mathematical model for the behavior of pedestrians[END_REF]; [START_REF] Helbing | On the fluid mechanics of human crowd motion[END_REF]) where pedestrian acts like automata (even with a memory criterium). Here this can be typically understood by looking at this forward-backward structure and the term in the divergence of the Kolmogorov equation ( 5), which is not fixed a priori since it is the feedback control.

• Moreover, the non-cooperative equilibrium approach is a bit restrictive since it does not allow to model more realistic pedestrian crowds with small social groups therein (2 or 3 people groups for instance).

• Finally a satisfying point we want to mention is that contrary to other game models, it allows entries and exits of pedestrians in the crowd (see section 4 for an example). Indeed, in games with a finite number N of players, it is very complicated to consider a nonconstant N.

In the next sections, we illustrate how MFG can be used for modeling pedestrian crowds. More precisely we introduce two different models for 2-population interactions, the first one is based on crowd aversion while the second one takes into account congestion in the running cost. We present numerical discretizations and simulations in both cases.

Pedestrian dynamics in the case of xenophobia: symmetry breaking

In this part we model two groups dynamics, each of whom consisting of infinitely many pedestrians. Consequently we use the finite horizon setting. We will focus here on a MFG that can be written as an optimal control problem (the running cost has the form f (t, x, α, m) = L(t, x, α) + V(x, m(x)) as mentioned in section 2. The optimal control formulation serves as a starting point for our numerical discretization, namely a gradient descent method. Note that we will only use the macroscopic optimal control formulation, but nevertheless, it is good to keep in mind the microscopic basis. The basic modeling assumption is that pedestrians in both groups have an aversion towards crowds. Furthermore the groups want to avoid each other (xenophobia). This approach differs from the one we use in the next section where the pedestrians' choices take congestion into account (that is the cost to move increases with the crowd distribution). Here crowd aversion only consists in penalizing congestion, not modeling it. The first of the previous aspects is mathematically expressed as convexity in the group problem while we will see that the latter leads to the loss of convexity in the joint problem for both groups.

We start with the macro-Nash problem between two populations (each of them being MFG), then state some theoretical results, detail the numerical algorithm and present some simulations including symmetry breaking.

Writing a macro-Nash problem between two groups

Let us focus on the case where two populations interact inside a subset Ω of R 2 . We also define the time-space domain Q = [0, T ] × Ω. We want to study equilibria between the two groups (typically Nash points as suggested by [START_REF] Lasry | Mean field games[END_REF]). Formally, the global optimization problem (linked to a continuum of individual problems) of group i, given the control and the mass evolution of the other group (i.e. (α -i t , m -i t )), reads as:

inf α i J i λ (α)
where

J i λ (α) = Q |α i | 2 2 m i + T 0 Φ i λ (m 1 t , m 2 t ) + Ω Ψ i m i (T ). (11)
Note that the upper index now refers to the group (or population) and not to the individual player. From now on, m i depends upon α i , more precisely it is viewed as a bounded nonnegative measure (i.e. belonging to the set M b (Q, R + )) which is a weak solution of the Kolmogorov equation:

∂ t m i - σ 2 2 ∆m i = -div(α i m i ), m i (0) = m i 0 . (12)
We distinguish the populations by considering different initial densities m i 0 (.) and different final incentive costs Ψ i , for i = 1, 2. However we study the simple case where the Brownian motion and the noise are similar for both groups. In the definition of the criterion (11), the coupling cost we consider is

Φ i λ (m t ) := Ω (m i t ) 2 + λm 1 t m 2 t ,
for a nonnegative real constant λ. This models situations in which individuals have an aversion towards members of his/her own group as well as the other group. Therefore λ can be interpreted as the "xenophobia" parameter. To see that, note that the individual mean field criterion is:

V i [m t ](x) = 2m i t (x) + λm -i t (x).
In fact we consider that xenophobia is significant when the aversion λ to the other group is higher than the one to the own group (equal to 2). The Nash problem between the groups is then:

(N) Find α = (α 1 , α 2 ) such that:

J i λ (α) = inf α i ∈M b (Q,R d ) J i λ (α i , α -i ), for i = 1, 2.

Existence and optimality 3.2.1. Optimality conditions

In this part we present a characterization of the Nash equilibria. To do so, let us introduce the MFG system for two groups: for i = 1, 2,

∂ t m i - σ 2 2 ∆m i + div(m i ∇v i ) = 0, m i (0) = m i 0 , (13)
∂ t v i + σ 2 2 ∆v i + |∇v i | 2 2 = Φ i λ (m) ′ , v i (T ) = Ψ i , (14)
and the joint minimization problem (for the two groups)

(Q) inf α=(α 1 ,α 2 ) J λ (α) := J 1 λ/2 (α) + J 2 λ/2 (α),
under the constraints: m i is a solution of ( 12), for i = 1, 2. Note that J λ is convex if the xenophobia parameter is small, that is λ ≤ 2. Now we are able to state the following optimality conditions.

Proposition 3.1. If λ ≤ 2 then the following assertions are equivalent:

1. α ∈ M b (Q, R d ) is a solution of (N) and m satisfies (12) for α = α, 2. α ∈ M b (Q, R d
) is a solution of (Q) and m satisfies ( 12) for α = α, 3. (m, v) is a solution of the MFG system ( 13)-( 14), with α = α = ∇v.

If λ > 2 then it is only necessary i.e. 2. ⇒ 1.,3.

A proof of this statement is based on classical arguments of differential calculus and can be found in [START_REF] Lachapelle | Control and transport problems in economics: Theoretical and numerical aspects[END_REF]. It basically consists of introducing first order perturbations, the linearized PDE and compute the first variation.

Existence

If problem (Q) has a solution, then proposition 3.1 ensures the existence of a Nash point between the two groups. In order to give an existence result, we may reformulate the problem following [START_REF] Benamou | A computational fluid mechanics solution to the monge-kantorovich mass transfer problem[END_REF]; [START_REF] Buttazzo | An optimization problem for mass transportation with congested dynamics[END_REF]; [START_REF] Lachapelle | Computation of mean field equilibria in economics[END_REF]. We adopt a vectorial point of view, and use the following notations:

• m = (m 1 , m 2 ), and for all x = (x 1 , x 2 ) ∈ R 2 + * , 1 x := ( 1 x 1 , 1
x 2 ), • q = (q 1 , q 2 ), and for all y = (y 1 , y

2 ) ∈ R d × R d , |y| 2 := (|y 1 | 2 , |y 2 | 2 ), • A := {(q, m) : Q (∂ t u+ σ 2 2 ∆u)dm i + Q ∇u.dq i = Ω Ψ i m i (T )-u 0 m i 0 , ∀u ∈ C ∞ (Q), i = 1, 2}.
Next we formally introduce the new variables q = (αm, m) and the corresponding cost functional

K(q, m) = J λ (α), if q << m and q = αm +∞ else,
where q << m means that the measure q is absolutely continuous with respect to m. In this setting, a rigorous formulation of (

Q) is inf (q,m)∈A K(q, m). ( 15
)
The reformulation allows us to prove the following proposition.

Proposition 3.2. If λ ≤ 2 and m 1 0 , m 2 0 ∈ L 2 (Q), then problem (Q) possesses a solution (which is unique as soon as λ < 2). Moreover there exists a Nash point i.e. a solution of (N).

Remark 5. Note that existence does not fail when adding a constraint of the type m ≤ constant as in Maury et al. (to appear).

Given the reformulation ( 15), the proof is a simple adaptation of the one obtained by [START_REF] Buttazzo | An optimization problem for mass transportation with congested dynamics[END_REF]. The key points are to work on the space of bounded measures, to prove that the functional K is lower-semicontinuous for the weak-⋆ topology and use Jensen's inequality. A complete proof of this proposition is provided in [START_REF] Lachapelle | Control and transport problems in economics: Theoretical and numerical aspects[END_REF]. In the next section we deal with defining a numerical procedure to approximate the solution.

Numerical setting

In this part we introduce the discretization and a gradient descent method in order to approximate the solution(s) of problem (N). More precisely, we distinguish the cases when the joint problem (Q) is convex from when it is not. In the convex setting (i.e. when λ ≤ 2) we describe the gradient descent that we apply to the joint functional. The non-convex case λ > 2 (in which the xenophobia is significant) is more involved but interesting (we expect non-uniqueness). Here we use an alternating directions method which takes advantage of the convexity of group i's problem, given group (-i)'s evolution.

Gradient First of all, let us write the gradient formula of the functional. We look at the reformulated problem given by ( 15). We slightly modify the point of view considering that the density m is an affine function of the momentum q. To fix ideas, the joint problem reads as:

inf q F(q) := i=1,2 Q |q i t | 2 2m i t + Φ i λ/2 (m t ) + Ω Ψ i m i T , (16)
where m i , i = 1, 2 solves :

∂ t m i - σ 2 2 ∆m i = -div(q i ), m i (0, .) = m i 0 (.). ( 17
)
The gradient of the functional can be calculated using classical differential calculus,

∀(q, m) ∈ A, ∀w = (w 1 , w 2 ) ∈ M b (R 2d ), ∇F(q).w = Q q i m i + ∇θ i .dw i i=1,2 , (18)
where θ i satisfies, for i = 1, 2:

θ i | t=T = Ψ i and -∂ t θ i - σ 2 2 ∆θ i = - |q i | 2 2(m i ) 2 + (2m i + λm -i ). ( 19
)
Algorithm for the convex case Since problem ( 16) is convex if λ ≤ 2, we use a gradient descent method. We focus on the 2D-case (d = 2) and take Ω = [0, 1] 2 with periodic boundary conditions. Let M and N be two positive integers, we define the time and space steps by dt = 1 N and dx = 1 M . For (i, j, k) ∈ A := {0, ..., N} × {0, ..., M} 2 , for a given function f defined on Q, f i j,k denotes the numerical approximation of f (idt, jdx, kdy). Equations ( 17) and ( 19) are iteratively solved by a finite difference method. We use the following approximations

∂ t f (idt, jdx, kdy) = f i+1 j,k -f i j,k dt , ∆ f (idt, jdx, kdy) = f i j+1,k -2 f i j,k + f i j-1,k (dx) 2 + f i j,k+1 -2 f i j,k + f i j,k-1 (dy) 2 . At step n, let f (n) := f i,(n) j,k (i, j,k)∈A
. Then the gradient descent method (GDM) reads as:

1. Initialization: Choose q (0) then compute m (0) by solving (17) with the finite difference scheme. 2.

Step n:

• Compute θ (n) by solving numerically (19) with q (n-1) and m (n-1) , then compute the discretized gradient ∇F(q (n-1) ) (formula (18)), using θ (n) . • Compute the descent: q (n) = q (n-1)ρ n ∇F(q (n-1) ).

• If ||q (n)q (n-1) || <Tol1, then stop the algorithm (Tol1 is a tolerance threshold defined by the user).

Else, n = n + 1.

Note that ρ n above is the descent step size, it is chosen optimal, i.e. minimizing the following: ρ ∈ [0, 1] → F(q (n-1) -ρ∇F(q (n-1) )).

Alternating directions method for the non-convex problem We expect non-uniqueness in the case where aversion to the other group is significant (λ > 2). It is well known that classical gradient descent methods have convergence problems or diverge due to the non-convexity of the functional. Therefore we use the alternating directions method, which uses the convexity of both group i's problem. The algorithm reads as: Given group (-i)'s evolution solve

inf q i F i (q) := T 0 Ω |q i t | 2 m i t + Φ i λ (m t )dt + Ω Ψ i m i T ,
where m i solves (17) for q i , i = 1, 2. One can easily get the formula of the gradient of F i looking at the joint case ( 18)-( 19). Note that we can deduce m from our computation of q. The idea of the alternating directions method is to apply GDM successively for each group. Note that the upper index refers to the group number and the lower one to the iteration.

Initialization:

Choose q 1 0 then compute q 2 1 with GDM and q 1 0 .

2.

Step k 1: We know q 2 k . • Compute q 1 k then q 2 k+1 by using GDM (with, respectively, q 2 k and q 1 k ).

• If ||q i k -q i k-1 || <Tol2 for i = 1, 2, then stop the procedure. Else, k = k + 1.

Simulations

The GDM shows good convergence results when the initials densities of individuals are significantly positive (i.e. m i 0 > constant > 0). We set T = 1 and σ 2 2 = 0.01 in all numerical simulations in this section.

Test 1: crowd aversion in a single group

In the first example we focus on a case involving only one population (m 2 0 = 0), i.e. a similar framework as the one studied by [START_REF] Buttazzo | An optimization problem for mass transportation with congested dynamics[END_REF]. Fig. 1 shows the initial density of agents (centralized around the point (0.1, 0.1)) and the final cost, modeling a strong incentive for individuals to be in some neighborhood of (0.5, 0.8) and (0.8, 0.5) at instant T . Note that the running time of dispersion is greater than one half. We then observe in Fig. 3 a split inside the population so that individuals can converge to the two attractive areas. Finally, the discrete energy seems to reach quickly the minimum (5 iterations), see Fig. 4.

Test 2: groups interactions

Next we consider a more interesting example with two populations. Recall that we look for Nash equilibria between two groups whose global optimization problem is (15). We use the procedure detailed before (starting with group 2). In order to emphasize the xenophobia behavior we choose λ = 20 in the definition Φ i λ (m 1 , m 2 ) = Ω (m i ) 2 + λm 1 m 2 . We consider a symmetric configuration and represent the graphs of the initial densities (m i 0 , i = 1, 2) and final costs (Ψ i , i = 1, 2) in Fig. 5. Group 1 is initially centered around (0.35, 0.5), group 2 around (0.5, 0.35). Concerning the final costs the situation is still symmetric since they model incentives to reach positions in the neighborhood of (respectively for group 1 and 2) (0.65, 0.5) and (0.65, 0.5). With such a situation we are interested in crossing phenomenon. The graphs of both densities are depicted in Fig. 6. Note that we observe the same spreading as in Test 1. However, the most interesting evolution period is described in Fig. 7. Indeed, we can see that group 1 gives the priority to group 2 to go to its attractive area passing through the center of the domain (the shortest road for the euclidian metric). Some of the individuals of group 1 wait, some others go through the border (periodic conditions), and the lasts go through the center (the most congested area). Anyway we note that group 2 reaches its goal quicker than group 2. Looking at Fig. 8, we can check that both group are finally centered around the points (0.65, 0.5) and (0.5, 0.65). Finally we remark that we obtain the opposite (or symmetric) situation, i.e. when starting to optimize with group 1. Then, the symmetry breaking seems to be a consequence of this starting choice. This section is devoted to mean field game congestion models. We start with the motivation for a single population, and then present the generalization for two populations. Note that although the basic assumptions on the cost functional are quite simple, the model is capable of reproducing interesting phenomena like lane formation.

Single species congestion model

Here we define as in Lions (2008Lions ( -2010) ) the running cost as f i (t, x, α, m) = ( |α| q q ( x-i ) a + k(t, x i ))e -rt , where:

• k models panic (this permanent cost does not depend on density or time, the larger it is the faster people want to get out),

• r > 0 is the discount factor (i.e. the actual position of all other people affects the strategy of an individual much more than the evolution in time),

• we write here the general case where q > 1 but we will in fact focus, as usual, on the quadratic case q = 2

• the nonnegative power index a models how congestion impacts the cost to move. Note that the running cost increases with the a-th power of the distribution of people, forcing individuals to avoid congested areas. Note that for a = 0 the cost does not penalize congestion (for a fixed value of q). Furthermore we can show uniqueness of solutions if a < a * , where a * is a nonnegative constant depending on q (see Lions (2008Lions ( -2010))) for more details).

Then the corresponding stationary MFG (see [START_REF] Lasry | Mean field games[END_REF] and Lions (2008Lions ( -2010) ) for the detailed derivation) is given by

ν∆u + 1 p |∇u| p m b + -ru + λ = k (20a) -ν∆m + div(m (∇u) p-1 m b ) = 0 (20b) m > 0, udx = 0, m dx = 1, (20c)
with ν = σ 2 2 , b = a q-1 and p = q q-1 .

Two-population model

The corresponding two species model looks slightly different. Here we assume that each density would like to avoid congestion within its own group as well as with the other. The corresponding running cost is given by

f i (t, x, α, m) = |α| q q (m i ) a (m j ) ã + k(x, t) e -rt (21)
for i = 1, 2, with m = (m 1 , m 2 ). Note that when ã > a, individuals of one crowd primarily avoid congestion with people from the other group. Then the corresponding mean field game for m 1 and m 2 reads as

ν∆u i + 1 p |∇u i | p m b i m b j -ru i + λ i = k (22a) -ν∆m i + div(m i (∇u i ) p-1 m b i m b j ) = 0 (22b) m i > 0, u i dx = 0, m i dx = 1, (22c)
with i, j = 1, 2 and b = ã q-1 . To avoid the division by zero we consider a slightly different model (which corresponds to consider (c + m) instead of m in the running cost ( 21)), namely

-ν∆u i + 1 p |∇u i | p (c + m i ) b (c + m j ) b -κu i + λ i = k (23a) -ν∆m i -div(m i (∇u i ) p-1 (c + m i ) b (c + m j ) b ) = 0 (23b) m i > 0, u i dx = 0, m i dx = 1. ( 23c
)
for a positive constant c. We choose homogeneous Dirichlet conditions for m i at the exit (people leave the room, hence the density has to be zero) and a homogenous Neumann boundary conditions on the rest of the boundary, i.e. m i = 0 for all x ∈ Γ out i and ∂m i ∂n = 0 on the rest of the boundary.

The same boundary conditions are set for u i . Here the additional variable λ i and the integral condition for u i are not necessary. We replace the integral condition for m i by a source term in the Kolmogorov equation, i.e.

-ν∆m i -div(m i (∇u i ) p-1 (c + m i ) b (c + m j ) b ) = f (x).
This source term can be interpreted as an exit of an underground or a large supermarket.

Figure 9: Sketch of the computational domain with two exits

Numerical simulations

Finally we would like to illustrate the behavior of the proposed congestion model ( 23) with various numerical simulations. Note that system ( 23) is a coupled system of four nonlinear PDEs, which can be solved using Newtons' method. The abstract formulation of ( 23) is given by

F(u 1 , m 1 , u 2 , m 2 ) = 0
where F is a nonlinear operator. Hence Newtons' method reads as

JF(w n )q n = -F(w n ), w n+1 = w n + τq n , where w n = (u n 1 , m n 1 , u n 2 , m n
2 ) denotes the n-th Newton iterate, JF(w n ) is the Jacobi matrix evaluated at w n and τ is an appropriate damping parameter. We use a hybrid discontinuous Galerkin method to discretize the Newton system in space. All computations and simulations are based on the softwares package Netgen/NgSolve (cf. [START_REF] Schöberl | Netgen-an advancing front 2d/3d-mesh generator based on abstract rules[END_REF]) and Pardiso (cf. Schenk and Gärtner (2004b,a)). Throughout this section we choose the following parameters q = 2, r = 1 and k = 1.

Validation of the model

First we would like validate our mathematical model with a simple test case. We choose a rectangular domain of size [-1, 1] × [-0.2, 0.2], that has two exits (see Figure 9). We have two sources f i , i = 1, 2 (Gaussians) located at (0.8, -0.1) and (-0.8, -0.1). The first group located at (0.8, 0.1) wants to exit at the upper left door, while the other one wants to get to the lower right one. We assume that people want to avoid congestion within their own and the other group at the same level, i.e. a = ã = 0.5. The diffusivity parameter is set to ν = 0.05. The intuitive solution of this problem is that each group stays in its "predefined" lane, that is determined by the offset of both sources. Our numerical simulations confirm this intuitive solution, we observe the formation of two lanes in Figure 10.

Avoidance behavior

Next we would like to study the avoidance behavior if two groups move towards each other without the possibility to totally avoid the other group. We choose the same domain as in the first example, only the exits are different. Here the exits are located at x = ±1. We put two sources of people at (±0.8, 0), each group wants to move to the opposite exit. Again we set a = ã = 0.5 and ν = 0.05. The solution is depicted in Figure 11, we observe that both groups move towards each other. Because of the high densities at x = ±0.8, the groups split and merge behind the high density regions.

Obstacle problem -symmetry breaking

Let us consider a domain with a circular obstacle, i.e. a rectangle [-1.5, 1.5] × [-0.2, 0.2] with a circular obstacle of radius r = 0.05 located at (0, 0). Here we place two sources of people at (±0.75, 0), the exists are located on the right and the left of the domain. We assume that the groups do not strongly mind congestion within themselves, but clearly with the other group. Hence we choose two different coefficients a and ã, i.e. a = 0.25 and ã = 0.75. Here the diffusion parameter is set to ν = 0.05. The solution of the Newton iteration is illustrated in Figure 12. Here we observe the loss of symmetry and the consequent formation of lanes, i.e. each group passes the obstacle on one side. We think that the different values of a and ã facilitate the symmetry breaking. Furthermore we suspect that the geometry of the domain as well as the diffusion parameter play another role in this complex process.

Corridors

Finally we consider a more complicated domain, i.e. a rectangular domain which has an L-shaped obstacle in the center (see Figure 13) . This obstacle creates a small "door" on the left side. We place two sources at the lower left and lower right corner at (0.35, 0.3) and (1.4, 0.3). The exits are located at the right and left side of the rectangle, each group wants to get to the opposite exit. We set a = 0.5 and ã = 2, i.e. people within a group try to avoid congestion with the other group but do not mind congestion as much within their own group. The diffusion parameter is set to ν = 0.1. Again we observe lane formation, the first group takes the small door to exit at the top edge, while the second group takes the left corridor to exit on the bottom edge (see Figure 13). From our viewpoint this is a reference example where pure forward models may fail to describe this "anticipation" solution. Indeed, there is no reason why the pedestrians coming from the right source should avoid the lower corridor if they do not anticipate that people coming from the left source will reach the small door much faster (which is naturally a congested area).

Conclusion and outlook

In this paper we presented a new mathematical modeling approach for crowd motion models. It is based on the theory of mean field games, which assumes that individuals are smart players which try to optimize their strategy and path with respect to certain costs (equilibrium with rational expectations approach). This clearly differentiates MFG from robot or automata models. We focused on modeling aversion and congestion behavior between two groups of pedestrians via MFG. In both cases we provided numerical solutions to approximate the solutions. Here we were able to reproduce characteristic features of crowds like lane formation or symmetry breaking. We showed that many situations can be treated via MFG (finite horizon vs stationary problems, aversion vs congestion, objective function vs panic, boundary conditions and sources terms). The main feature of MFG is the forward-backward structure, the usual forward part describes the crowd evolution while the backward gives the process of how the expectations are built.

To the best of our knowledge this is one of the first applications of MFG to pedestrian crowd modeling. It is thus necessarily imperfect and many things remain to be done. Here we list some of them. First it would be interesting to introduce a common noise (e.g. modeling an earthquake incidence or a bridge failure's rumor etc.). This is feasible in a MFG setting, but leads to infinite dimensional PDEs (thus more complicated). Furthermore we would like to include the realistic assumption that people can not observe the whole distribution of the crowds. This could be accomplished by introducing partial blindness. The goal is then to model myopic individuals. We think that the proposed MFG models give a good impression of the possibilities of mean field 20 games theory. These models may serve as a first and promising step towards the development of more realistic MFG models for crowd motion.

Figure 2 :

 2 Figure 1: Data

  Figure 3: Splitting and centralization during the second half

  Final costs Ψ i

FigureFigure 6 :Figure 7 :

 67 Figure 5: Data

Figure 11 :

 11 Figure 10: Validation: top view of both groups

Figure 12 :

 12 Figure 12: Obstacle problem: top view of both groups

Acknowledgements

Both authors would like to thank Jean-Michel Lasry for his helpful comments and suggestions. Furthermore we would like to thank Gabriel Turinici and Peter A. Markowich for initiating this work. MTW acknowledges financial support by the Austrian Science Foundation FWF via the Hertha Firnberg Project T456-N23.


Emmanuel Fricain
email: fricain@math.univ-lyon1.fr

Javad Mashreghi
email: javad.mashreghi@mat.ulaval.ca

ON A CHARACTERIZATION OF FINITE BLASCHKE PRODUCTS

Keywords: 2000 Mathematics Subject Classification. Primary: 30D50, Secondary: 32A70 Blaschke products, zero sets, convergence, automorphism

We study the convergence of a sequence of finite Blaschke products of a fix order toward a rotation. This would enable us to get a better picture of a characterization theorem for finite Blaschke products.

Introduction

Let (z k ) 1≤k≤n be a finite sequence in the open unit disc D. Then the rational function

B(z) = γ n k=1 z k -z 1 -zk z ,
where γ is a unimodular constant, is called a finite Blaschke product of order n for D [START_REF] Mashreghi | Representation Theorems in Hardy Spaces[END_REF]. There are various results characterizing these functions. For example, one of the oldest ones is due to Fatou.

Theorem A (Fatou [START_REF] Fatou | Sur les fonctions holomorphes et bornées à l'intérieur d'un cercle[END_REF]). Let f be analytic in the open unit disc D and suppose that

lim |z|→1 |f (z)| = 1.
Then f is a finite Blaschke product.

For an analytic function f : Ω 1 -→ Ω 2 , the number of solutions of the equation

f (z) = w, (z ∈ Ω 1 , w ∈ Ω 2 ),
counting multiplicities, is called the valence of f at w and is denoted by v f (w). It is well-known that a finite Blaschke product of order n has the constant valence n for each w ∈ D. But, this property in fact characterizes finite Blaschke products of order n.

Theorem B (Fatou [START_REF] Fatou | Sur les équations fonctionnelles[END_REF][START_REF] Fatou | Sur les équations fonctionnelles[END_REF][START_REF] Fatou | Sur les équations fonctionnelles[END_REF], Radó [9]). Let f : D -→ D be an analytic function of constant valence n ≥ 1. Then f is a finite Blaschke products of order n.

Our main concern in this paper is the following result of Heins. We remind that a conformal automorphism of D has the general form

T a,γ (z) = γ a -z 1 -ā z ,
where a ∈ D and γ ∈ T. Instead of T a,1 we simply write T a . The collection of all such elements is denoted by Aut(D). The special automorphism

ρ γ : D -→ D z -→ γz is called a rotation. Note that ρ γ = T 0,-γ .
Let us remind two further notions. In the following, when we say that a sequence of functions on D is convergent, we mean that it converges uniformly on compact subsets of D.

If f : D -→ C and (a k ) k≥1 ⊂ D, with lim k→∞ |a k | = 1, is such that lim k→∞ f (a k ) = L, then L is called an asymptotic value of f .
Theorem C (Heins [START_REF] Heins | Some characterizations of finite Blaschke products of positive degree[END_REF]). Let f : D -→ D be an analytic function. Then the following assertions are equivalent:

(i) f is a finite Blaschke product of order ≥ 1;

(ii) if the sequence of automorphisms T a k ,γ k , k ≥ 1, tends to a constant of modulus one, and if

g k = T f (a k γ k ) • f • T a k ,γ k is convergent, then g k tends to a rotation; (iii) we have lim |z|→1 (1 -|z| 2 ) |f ′ (z)| 1 -|f (z)| 2 = 1;
(iv) f has no asymptotic values in D and the set {z ∈ D : f ′ (z) = 0} is finite.

Our contribution is to further clarify the item (ii) above in which it is assumed that the limits of two sequences of analytic functions exist. No doubt the existence of these limits depend on the parameters a k ∈ D and γ k ∈ T, k ≥ 1.

It is easy to see that (T a k ,γ k ) k≥1 tends to a unimodular constant γ 0 if and only if

(1.1) lim k→∞ a k γ k = γ 0 .
In the first place, we show that if B is a finite Blaschke product of order ≥ 1 and (1.1) holds, then

T B(a k γ k ), γk • B • T a k ,γ k (z) -→ B ′ (γ 0 ) |B ′ (γ 0 )| z
as k -→ +∞ and the convergence is uniform on compact subsets of D. See Theorem 3.2. Therefore, the sequence T B(a k γ k ) • B • T a k ,γ k , which was considered in Theorem C, is convergent to a rotation if and only if γ k is convergent. We refer to the example given at the end to see that this condition cannot be relaxed. In this special situation, if a k -→ γ 0 , as k -→ +∞, then

T B(a k ) • B • T a k (z) -→ B ′ (γ 0 ) |B ′ (γ 0 )| z,
uniformly on compact subsets of D. The above observations enable us to rewrite Theorem C in the following form.

Theorem D. Let f : D -→ D be an analytic function. Then the following assertions are equivalent:

(i) f is a finite Blaschke product of order ≥ 1; (ii) if γ k a k -→ γ 0 ∈ T, then T f (a k γ k ),γ k • f • T a k , γ k tends to a rotation; (iii) if a k -→ γ 0 ∈ T, then T f (a k ) • f • T a k tends to a rotation; (iv) the equality lim |z|→1 (1 -|z| 2 ) |f ′ (z)| 1 -|f (z)| 2 = 1 holds; (v) f has no asymptotic values in D and the set {z ∈ D : f ′ (z) = 0} is finite; (vi) f has a constant valence on D.
Moreover, if the above conditions hold, then the rotation promised in (ii) and (iii)

is ρ γ , where γ = f ′ (γ 0 )/|f ′ (γ 0 )|.
We will just study the implication (i) =⇒ (ii) in Theorem 3.2. That (ii) =⇒ (iii) is trivial. The implication (iii) =⇒ (iv) is an easy consequence of the formula

(1 -|a k | 2 ) f ′ (a k ) 1 -|f (a k )| 2 = T f (a k ) • f • T a k ′ (0),
and that, by assumption, the latter tends to a unimodular constant. The more delicate steps (iv) =⇒ (v) =⇒ (vi) and (vi) =⇒ (i) are already taken respectively by Heins and Fatou-Radó.

In the course of proof, we also show that, for a finite Blaschke product B, the zeros of B ′ which are inside D are in the hyperbolic convex hull of the zeros of B.

hyperbolic convex hull

The relation

B(z) B(1/z) = 1 shows that B ′ (z) = 0 ⇐⇒ B ′ (1/z) = 0.
Hence, to study the singular points of B, it is enough to consider the zeros of B ′ which are inside the unit disc D. Note that B has no singular points on T, according to (3.1) for instance. Moreover if B has n zero in D, then it is easy to see that B ′ has n -1 zeros in D (counting with multiplicities).

There are various results about the relations between the zeros of a polynomial P and the zeros of its derivatives. The eldest goes back to Gauss-Lucas [START_REF] Lucas | Propriétés géométriques des fractionnes rationnelles[END_REF] which says that the zeros of P ′ are in the convex hull of the zeros of P . Recently, Cassier-Chalendar [START_REF] Cassier | The group of the invariants of a finite Blaschke product[END_REF] established a similar result for Blaschke products: the zeros of B ′ in D are in the convex hull of the zeros of B and {0}. We show that the zeros of B ′ are in the hyperbolic convex hull of the zeros of B. We do not need to enlarge the zero sets of B by adding {0}. Moreover, the hyperbolic convex hull of a set is a subset of the Euclidean convex hull of the set and {0}. In particular, we improve the result obtained in [START_REF] Cassier | The group of the invariants of a finite Blaschke product[END_REF].

Let z 1 , z 2 ∈ D. Then the hyperbolic line between z 1 and z 2 is given by

[0, 1] -→ D t -→ z 1 -z 1 -z 2 1-z 1 z 2 t 1 -z 1 z 1 -z 2 1-z 1 z 2 t .
This representation immediately implies that the three distinct points z 1 , z 2 , z 3 ∈ D are on the same hyperbolic line if and only if

(2.1) z 1 -z 2 1 -z 1 z 2 / z 1 -z 3 1 -z 1 z 3 ∈ R.
Furthermore, we see that the hyperbolic lines in D can be parameterized by

γ a -z 1 -az = t, t ∈ [-1, 1],
with γ ∈ T and a ∈ D.

Adopting the classical definition from the Euclidean geometry, we say that a set

A ⊂ D is hyperbolically convex if z 1 , z 2 ∈ A =⇒ ∀t ∈ [0, 1], z 1 -z 1 -z 2 1-z 1 z 2 t 1 -z 1 z 1 -z 2 1-z 1 z 2 t ∈ A.
The hyperbolic convex hull of a set A is the smallest hyperbolic convex set which contains A. This is clearly the intersection of all hyperbolic convex sets that contain A. Of course, if we need the closed hyperbolic convex hull, we must consider the intersection of all closed hyperbolic convex sets that contain A.

In the proof of the following theorem, we need at least two elementary properties of the automorphism Proof. First suppose that all zeros of B are in D + . Then, taking the logarithmic derivative of B, we obtain

T a (z) = a -z 1 -ā z , (a, z ∈ D).
B ′ (z) B(z) = n k=1 1 -|z k | 2 (1 -zk z)(z -z k ) , which implies ℑ B ′ (z) B(z) = n k=1 ℑ 1 -|z k | 2 (1 -zk z)(z -z k ) .
In the light of last expression, put

Φ(z) = 1 -|a| 2 (1 -ā z)(z -a) ,
where a D + is fixed. We need to obtain the image of D -under Φ. To do so, we study the image of the boundary of ∂D -under Φ. On the lower semicircle

T -= {e iθ : -π ≤ θ ≤ 0} we have Φ(e iθ ) = 1 -|a| 2 (1 -ā e iθ )(e iθ -a) = 1 -|a| 2 |e iθ -a| 2 e -iθ .
Therefore, T -is mapped to an arc in C + . Moreover, on the interval t ∈ (-1, 1), we have

Φ(t) = 1 -|a| 2 (1 -ā t)(t -a) = 1 -|a| 2 |t -a| 2 × t - ā 1 -ā t = 1 -|a| 2 |t -a| 2 T t (ā).
Thus (-1, 1) is also mapped to an arc in C + = {z ∈ C : ℑz > 0}. In short, ∂D - is mapped to a closed arc in C + . Therefore, we deduce that Φ maps D -into C + . Note that 0 is not contained in Φ(D -). This fact implies that B ′ does not have any zeros in D -. By continuity, we can say that if all zeros of B are in D + , then so does all the zeros of B ′ which are in the disc.

A general automorphism of D has the form T a,γ = γT a , for some γ ∈ T and a ∈ D. Let f = B • T a,γ . Then f is also a finite Blaschke product with zeros T a (γ z 1 ), T a (γ z 2 ), . . . , T a (γ z n ). Moreover, if we denote the zeros of B ′ in D by w 1 , w 2 , . . . , w n-1 , the zeros of f ′ in D would be T a (γ w 1 ), T a (γ w 2 ), . . . , T a (γ w n-1 ).

Considering the boundary curves of the hyperbolic convex hull of z 1 , . . . , z n , if we choose a and γ such that

T a (γ z 1 ), T a (γ z 2 ), . . . , T a (γ z n ) ∈ D + ,
then the preceding observation shows that T a (γ w 1 ), T a (γ w 2 ), . . . , T a (γ w n-1 ) ∈ D + .

Therefore, we see that if the zeros of B are on one side of the hyperbolic line

γ a -z 1 -ā z = t, t ∈ [-1, 1],
then the zeros of B ′ are also on the same side. The intersection of all such lines gives the hyperbolic convex hull of the zeros of B.

Remark:

The argument above also works for infinite Blaschke products. Theorem 2.1 is sharp in the sense which is crystalized by the following examples. Let a, b ∈ D and put

B(z) = a -z 1 -ā z m b -z 1 -b z n .
Clearly B ′ has m + n -1 zeros in D which are a, m -1 times, and b, n -1 times, and the last one c is given by the equation

m(1 -|a| 2 ) (1 -ā z) 2 a -z 1 -ā z m-1 b -z 1 -b z n + a -z 1 -ā z m b -z 1 -b z n-1 n(1 -|b| 2 ) (1 -b z) 2 = 0.
Rewrite this equation as

z -a 1 -a z / z -b 1 -b z = - m(1 -|a| 2 ) |1 -ā z| 2 / n(1 -|b| 2 ) |1 -b z| 2 ,
which, by (2.1), reveals that a, b, c are on the same hyperbolic line. Moreover, as m and n are any arbitrary positive integers, the point c traverses a dense subset of the hyperbolic line segment between a and b.

With more delicate calculations, we can consider examples of the form

B(z) = a -z 1 -ā z m b -z 1 -b z n c -z 1 -c z p ,
where a, b, c ∈ D and m, n, p ≥ 1, and observe that the zeros of B ′ , for different values of m, n, p, form a dense subset of the hyperbolic convex hull of a, b, c.

Uniform convergence on compact sets

The following result is not stated in its general form. We are content with this version which is enough to obtain our main result that comes afterward.

Lemma 3.1. Let B(z) = γ n k=1 z k -z 1 -zk z and fix any M such that max{ |z 1 |, |z 2 |, . . . , |z n |} < M < 1.
Then there is a constant δ = δ(M, B) > 0 such that, for any two distinct points z, w in the annulus {z : M ≤ |z| ≤ 1/M}, we have

|z -w| < δ =⇒ B(z) = B(w).
Proof. Suppose this is not true. Then, for each k ≥ 1, there are a k and b k such that

M ≤ |a k |, |b k | ≤ 1/M, a k = b k , |a k -b k | < 1/k, but B(a k ) = B(b k ).
Since the closed annulus is compact, (a k ) k≥1 has a convergent subsequence. Without loss of generality, we may assume that a k -→ a for some a with M ≤ |a| ≤ 1/M. This implies b k -→ a. Therefore, we would have

B ′ (a) = lim k→∞ B(a k ) -B(b k ) a k -b k = 0.
This is a contradiction. as a matter of fact, by Theorem 2.1, the zeros of B ′ in D are in the hyperbolic convex hull of the zeros of B, and those outside D are the reflection of the zeros inside with respect to T. In other words, B ′ has no zeros on the annulus.

There is another way to arrive at a weaker version of Lemma 3.1, which will also suffice for us. On T, the nice formula

(3.1) |B ′ (e iθ )| = n k=1 1 -|z k | 2 |e iθ -z k | 2
holds. Hence, at least, B ′ (e iθ ) = 0 for all e iθ ∈ T. Thus, by continuity, there is an annulus on which B ′ has no zeros, and as we explained above, for such a region δ exists.

For the sake of completeness, we remind that if a k ∈ D, and γ k , γ 0 ∈ T, are such that a k γ k -→ γ 0 , as k -→ +∞, then the simple estimation

γ 0 -γ k a k -z 1 -āk z ≤ 2 |γ 0 -a k γ k | 1 -|z| , as k -→ +∞, reveals that (3.2) T a k ,γ k (z) -→ γ 0
uniformly on compact subsets of D.

Theorem 3.2. Let B be any finite Blaschke product of order ≥ 1. Let a k ∈ D and γ k , γ 0 ∈ T, be such that lim

k→∞ a k γ k = γ 0 . Then T B(a k γ k ), γk • B • T a k ,γ k -→ ρ γ ,
where γ = B ′ (γ 0 )/|B ′ (γ 0 )|. In other words,

T B(a k γ k ), γk • B • T a k ,γ k (z) -→ B ′ (γ 0 ) |B ′ (γ 0 )| z
as k -→ +∞ and the convergence is uniform on compact subsets of D.

Proof. Suppose that B has order n. Then, for each fixed k, the function

(3.3) f k = T B(a k γ k ), γk • B • T a k ,γ k
is also a finite Blaschke product of order n. Its zeros, say z k,1 , z k,2 , . . . , z k,n , are the solutions of the equation

B T a k ,γ k (z) = B(a k γ k ).
Certainly z = 0 is a solution. We need to show that the other zeros accumulate to some points on T. Denote the solutions of B(w) = B(a k γ k ) by

w k,1 , w k,2 , . . . , w k,n ,
and put w k,1 = a k γ k . As k -→ ∞, we have |w k,j | -→ 1, (1 ≤ j ≤ n),
and, by Lemma 3.1,

|w k,i -w k,j | ≥ δ, ( 1
≤ i = j ≤ n).
By the transformation z = T -1 a k ,γ k (w), the first solution is mapped to the origin. But the other solutions, thanks to the positive distance between them, are mapped to points close to T. More precisely, we have (3.4) lim

k→∞ |z k,j | = lim k→∞ |T -1 a k ,γ k (w k,j )| = 1, (2 ≤ j ≤ n). Write f k (z) = η k z n j=2 z k,j -z 1 -zk,j z ,
where η k ∈ T. This formula shows

f ′ k (0) = η k n j=2 z k,j and (3.3) reveals that f ′ k (0) = T ′ B(a k γ k ), γk (B(a k γ k )) B ′ (a k γ k ) T ′ a k ,γ k (0) = 1 -|a k γ k | 2 1 -|B(a k γ k )| 2 B ′ (a k γ k ).
Hence, we obtain the representation A simple example shows that in Theorem 3.2, T B(a k γ k ), γk cannot be replaced by T B(a k γ k ) . In other words, the rotation by γk is needed to obtain the convergence. To see this fact, let a k be any positive sequence on [0, 1) tending to 1. Let γ k = (-1) k , and put B(z) = z 2 . Then

f k (z) = 1 -|a k γ k | 2 1 -|B(a k γ k )| 2 B ′ (a k γ k )

Firstly,Theorem 2 . 1 .

 21 T a • T a = id. Secondly, if a ∈ (-1, 1), then T a maps D -= D ∩ {z : ℑz < 0} into D + = D ∩ {z : ℑz > 0}. In the same token, we define D + = D ∩ {z : ℑz ≥ 0}. Let B be a finite Blaschke product. Then the zeros of B ′ which are inside D are in the hyperbolic convex hull of the zeros of B.

-z 1 -Corollary 3 . 3 .

 133 zk,j z -→ 1, as k → ∞. It is also known that lim z→γ 0 1 -|B(z)| 2 1 -|z| 2 = |B ′ (γ 0 )|.We are done. Let B be a finite Blaschke product. Let a k ∈ D and γ 0 ∈ T be such that limk→∞ a k = γ 0 . Then T B(a k ) • B • T a k -→ ρ γ , (k -→ +∞),where γ = B ′ (γ 0 )/|B ′ (γ 0 )|.

T

  B(a 2k γ 2k ) • B • T a 2k ,γ 2k (z) -→ z while T B(a 2k+1 γ 2k+1 ) • B • T a 2k+1 ,γ 2k+1 (z) -→ -z. Hence, T B(a k γ k ) • B • T a k ,γ k is not convergent. Of course, either by Theorem 3.2 or by direct verification, T B(a k γ k ), γk • B • T a k ,γ k (z) -→ z,

This work was supported by NSERC (Canada), Jacques Cartier Center and ANR project FRAB (France).


Jérémie Botton
email: jeremie.botton@inserm.fr

Barbara Heude

Jean Maccario

Jean-Michel Borys

Agnes Lommez

Pierre Ducimetière

Marie Aline Charles

Parental body size and early weight and height growth velocities in their offspring

Keywords: growth, epidemiology, parental BMI and height

Background: Whereas weight or height at a given age are the results of the cumulative growth experience, growth velocities allows the study of factors affecting growth at given ages. Aim: to study the relationships between parental height and body mass index (BMI) and offspring's height and weight growths during infancy and childhood. Study design: From the FLVSII population-based study, 235 parent-child trios belonging to 162 families examined in 1999. Outcome measures: From medical records and previous FLVS examinations, child's height and weight history was reconstructed. Weight and height growth velocities from birth to seven years were estimated from a modelling of individual growth curve and correlated with parent's body size in 1999. Results: Ponderal index and length at birth was significantly associated with maternal but not paternal BMI and height. In the first six months, height growth velocity was significantly associated with maternal stature (at three months: 0.12±0.05 and 0.02±0.05 cm/mo for a 10 cm difference in maternal and paternal height respectively) and weight growth velocity with paternal BMI (at three months: 5.7±2.8 and 1.9±2.3 g/mo for a difference of 1kg/m 2 in paternal and maternal BMI respectively). Between two and five years, height growth velocity was more significantly associated with paternal height whereas weight growth velocity was more closely associated with maternal BMI. Conclusions: Early childhood growth is characterised by alternate periods associated specifically with maternal or paternal BMI and height. This novel finding should trigger the search for specific genetic, epigenetic or environmentally shared factors from the mothers and fathers.

Introduction

There is emerging evidence for periods during infancy and childhood during which children are more susceptible to develop fat mass [START_REF] Ekelund | Association of weight gain in infancy and early childhood with metabolic risk in young adults[END_REF][START_REF] Botton | Postnatal weight and height growth velocities at different ages between birth and 5 y and body composition in adolescent boys and girls[END_REF]. We recently reported results supporting the hypothesis of two critical windows in early childhood associated with an increased risk of obesity in adolescence: the first six months and from three to at least five years [START_REF] Botton | Postnatal weight and height growth velocities at different ages between birth and 5 y and body composition in adolescent boys and girls[END_REF]. From six months to three years, weight growth velocity was only associated with height in adolescence.

Determinants of weight and height growth specific to these periods need to be studied. A determinant consistently associated with children's weight and body mass index (BMI) is parental BMI. It is therefore of interest to assess how parent-child associations develop during different periods. Researchers have been interested for a long time in the relationships between parental characteristics and weight during the offspring's early life. In a study published in 1954, weight at birth and in early infancy were more correlated with maternal than paternal stature, but this difference disappeared before adulthood [START_REF] Cawley | Parental stature and birth weight[END_REF].

More recently, many studies have reported the relationships between maternal and sometimes paternal BMI with the offspring's BMI at birth or older ages: 3 years [START_REF] Whitaker | The relationship of childhood adiposity to parent body mass index and eating behavior[END_REF]; 7.5 years [START_REF] Smith | Is there an intrauterine influence on obesity? Evidence from parent child associations in the Avon Longitudinal Study of Parents and Children (ALSPAC)[END_REF]; 9 years [START_REF] Gale | Maternal size in pregnancy and body composition in children[END_REF]. They concluded that there was a positive association between mother's and offspring's BMIs; there were differing results about the relationship with the father's BMI.

Anthropometric relationships between parents and children throughout childhood have been reported in the FLVS study [START_REF] Heude | Anthropometric relationships between parents and children throughout childhood: the Fleurbaix-Laventie Ville Santé Study[END_REF]: BMI at birth was highly correlated with the mother's but not the father's BMI; correlations between offspring's and father's BMI were significant at two years.

To our knowledge, few studies have reported relationships between parental anthropometry and weight growth during infancy and childhood [START_REF] Blair | Family, socioeconomic and prenatal factors associated with failure to thrive in the Avon Longitudinal Study of Parents and Children (ALSPAC)[END_REF][START_REF] Griffiths | Differential parental weight and height contributions to offspring birthweight and weight gain in infancy[END_REF], and none considered both height and weight growth. The velocities allow to study factors affecting growth at given ages, whereas weight and height are the result of the prenatal and postnatal cumulative growth experience.

We aimed at studying the relationships between parents' BMI and height with both their offspring's weight and height growth velocities between birth and seven years.

Population and Methods

The subjects were recruited in Fleurbaix and Laventie, two neighbouring towns in northern France, with respectively 2488 and 4426 inhabitants in 1992, when the Fleurbaix-Laventie Ville Santé (FLVS) study started. The first part of the study (FLVS I), was a five-year followup of children participating in a nutritional education program at school and included 579 families [START_REF] Maillard | Macronutrient energy intake and adiposity in non obese prepubertal children aged 5-11 y (the Fleurbaix Laventie Ville Santé Study)[END_REF].

The second part (FLVS II) was an epidemiological study on the determinants of weight change in the population. This study was proposed to the families who participated in the FLVS I study. Written consent was obtained from children and parents. The study was approved by the ethics committee of Lille (France) and the computer files were declared to the "Commission Nationale Informatique et Liberté". We recruited and examined, between April and September 1999, 294 of the 393 families still living in these towns that could be contacted. All family members aged 8 years and over were examined at home in 1999 (1113 individuals). In this analysis, we included 235 parent-child trios belonging to 162 families with available data on parents' BMI and children anthropometry (weight and length or height at birth, nine months, two years and weight and height growth velocities). Other trios were not included because for 180 children the number of measurements during childhood failed to meet the requirement for growth modelling (see below). Among the other trios, paternal measurements were missing for 55, maternal for 4 and offspring's at birth, nine or two years for 41.

In offspring, weight and height from birth until twelve years were available from the "Carnet de santé", a booklet for parents and health practitioners recording the child's health data (precision of the measurements unknown), and as part of the longitudinal follow-up of the children in the FLVS I and II studies since the children were between five and twelve years (weight was measured in light clothes to the nearest 0.1 kg and height to the nearest cm). We calculated the ponderal index at birth as weight in kg over cubic length in meters, and BMI at nine months and two years as weight in kg over square height in meter. These indices are generally uncorrelated to length or height (confirmed in our data set).

In parents, anthropometric data were collected in 1999 by trained physicians during the FLVS II study. Weight was measured in light clothes to the nearest 0.1 kg, height to the nearest cm and BMI calculated.

Statistical analysis

We first estimated the weight and height growth velocities (a). Then, we used these velocities in models to test if there were differences between parents, and according to time (b).

a. Weight and height growth velocities estimation

The methodology used to estimate weight and height growth velocities has been described in detail elsewhere [START_REF] Botton | Postnatal weight and height growth velocities at different ages between birth and 5 y and body composition in adolescent boys and girls[END_REF]. Measurements of weight and height from birth to twelve years allowed us modelling each child's growth individually, using an adaptation of the Jenss growth curve model [START_REF] Berkey | Comparison of two longitudinal growth models for preschool children[END_REF] by adding a squared-time coefficient:

i i i i i i i ε + t) e + (d t c + t b + a = (t) y exp 2 
, where y i is the weight of the child i at age t, ε i being the residual. The same model was found to also fit the height growth.

The least-squares estimates were determined for each child using the SAS PROC NLIN. We included children with at least six measures available between birth and twelve years, two before two years and one after five years. The first derivative of the equation of the individual models calculated at selected age (0.25, 0.5, 1, 2, 3, 4, 5 and 7 years) gave the respective growth velocities. b. Descriptive and aetiological analyses General characteristics are given as means (SD) or geometric means [CI] for the parents' BMI. Differences between mother's and father's characteristics were tested using paired ttests. For children, all between gender comparisons used mixed effect models to take into account the family resemblance between brothers and sisters.

We estimated the relative contribution of each parent's height or BMI to the offspring's growth velocities by a multiple regression model at each selected age, including the parameter of both parents. The models took into account the family resemblance by introducing a random intercept and were additionally adjusted for mother's age, child's gender, and gestational age when parameter was measured at birth.. We also tested whether the differential effect of parental measures at these selected age (age-dependence) was significant, by testing the interaction term between parent's measures and age in a mixed model for longitudinal data: the variance-covariance structure which best fitted our data was an autoregressive heterogeneous one.

Differences in relationships between boys and girls were tested by including in the models an interaction term with offspring's gender. Finally, differences between father's and mother's contributions were tested by comparing the corresponding regression coefficients by a Wald test.

BMI being independent of height at the population level, we reported the relationships with parents' BMI rather than parents' weight.

Significance tests were two-sided and a p-value greater than 0.05 noted as non significant. We used SAS Software version 9.1.

Results

None of the anthropometric parameters of the 1999 FLVS II examination was significantly different in our study population between children included or not in the present analysis.

General characteristics (table 1)

In the 162 families, fathers had as expected higher weight, height and BMI than mothers. The 122 male offspring had a mean birth weight of 3.39 kg ± 0.53 versus 3.35 kg ± 0.52 for the 113 girls (p=0.14), and a mean birth length of 50.3 cm ± 2.3 versus 49.7 cm ± 2.3 in girls (p=0.003). The mean weight and height growth velocities at given ages during infancy and childhood are presented in figures 1A and 1B. Boys gained 0.67 ± 0.16 kg/mo at three months compared to 0.65 ± 0.13 kg/mo for girls (p=0.17). Three months height growth velocity was significantly greater in boys (p=0.02). After three months, there were no significant differences in height and weight growth velocities between boys and girls.

Parent's BMI: relationship with offspring's BMI and weight growth velocities (figure 2)

Parent's BMI was not significantly associated with either the offspring's height or height growth velocities (results not shown). Thus, the associations described above between parent's BMI and the offspring's weight growth are unlikely attributable to the offspring's height growth.

Offspring's BMI as endpoint. Mother's BMI was significantly associated with ponderal index at birth (figure 2A). A difference of 1 kg/m 2 in mother's and father's BMI was associated with a difference in ponderal index of 0.097 ± 0.04 kg/m 3 (p=0.01) and 0.046 ± 0.05 kg/m 3 (p=0.34) respectively. Both parental BMIs were associated with the offspring's BMI at nine months (difference of 0.055 ± 0.03 kg/m 2 and 0.069 ± 0.03 kg/m 2 in the child's BMI per kg/m 2 in the mother's (p=0.04) and father's BMI (p=0.04) respectively).

Offspring weight growth velocity as endpoint. Offspring weight growth velocity at three months was more strongly associated with the father's BMI (figure 2B: 5.7 ± 2.8 g/mo per kg/m 2 in father's BMI and 1.9 ± 2.3 g/mo per kg/m 2 in mother's BMI). Between six months and three years, parental BMIs were not associated with the offspring's weight growth velocity. At older ages, both parental contributions increased, more rapidly for the maternal contribution. These associations between father's or mother's BMI and offspring weight growth velocity were age-dependant (p-for-interaction with age=0.08 and 0.001 respectively).

Although interaction with offspring's gender was not significant, we observed that from three years onward, maternal contribution seemed preponderant in girls (supplementary figures). At seven years, the difference of the contributions between the two parents reached 2.9 ± 2.3 g/mo for 1 kg/m 2 (p=0.22) in girls and 0.13 ± 2.2 g/mo for 1 kg/m 2 (p=0.95) in boys.

Parental stature: relationship with offspring's height and height growth velocities (figure 3)

Offspring's length as endpoint. The offspring's birth length and height during infancy and childhood were significantly and positively associated with the mother's height over the whole period and with the father's height only at two years (figure 3A). At nine months, the test for a difference between maternal and paternal contribution was close to significance (p=0.08).

Offspring height growth velocity as endpoint. The mother's height was more strongly associated with the offspring's height growth velocity in the first six months of life than father's height (figure 3B: 0.12 ± 0.05 cm/mo at three months for a 10 cm difference in mother's height and 0.02 ± 0.05 cm/mo for 10 cm difference in father's height). At one year, the height growth velocity was not associated with parental stature. From two to five years, height growth was strongly associated with the father's height, whereas significant associations with the mother's height reappeared at four years and its effect remained much smaller than those with father's height (p-for-difference in contributions at three and four years=0.02 and 0.08 respectively). These associations between father's or mother's height and height growth velocity were age-dependant (p-for-interaction with age=0.13 and 0.02 respectively). In girls, at seven years the contribution of maternal height re-increases, becomes significantly higher than paternal contribution (supplementary figures) and significantly higher than in boys (p for interaction with offspring's gender =0.006).

Discussion

Height and ponderal index at birth were significantly associated with maternal but not paternal height and BMI. In the first six months of life, height growth velocity was associated with maternal stature and weight growth velocity with paternal BMI. Between two and five years, the reverse was observed: height growth velocity was more significantly associated with paternal height whereas weight growth velocity was more closely associated with maternal BMI.

Many authors highlighted the importance of parental determinants in the height and weight of their offspring [START_REF] Cawley | Parental stature and birth weight[END_REF][START_REF] Magnus | Paternal contribution to birth weight[END_REF][START_REF] Safer | Early body mass index and other anthropometric relationships between parents and children[END_REF][START_REF] Harvey | Parental determinants of neonatal body composition[END_REF]. Height and weight mirror multiple convergent genetic, environmental and hormonal interactions. However, it is currently admitted that height is mainly determined by genetic aspects in the absence of prolonged deficit in energy intake [START_REF] Palmert | Genetic approaches to stature, pubertal timing, and other complex traits[END_REF][START_REF] Silventoinen | Heritability of adult body height: a comparative study of twin cohorts in eight countries[END_REF]. Weight might be less heritable than height, even if genetic factors are probably still preponderant during at least infancy and childhood [START_REF] Van Dommelen | Genetic study of the height and weight process during infancy[END_REF]. The mother's contribution to the offspring's length at birth probably reflects both genetic factors and the intra-uterine constraint. After birth, it is thought that the height growth velocity adjusts towards a genetically predicted trajectory [START_REF] Veldhuis | Endocrine control of body composition in infancy, childhood, and puberty[END_REF]. However, our data show that the parental influences on the offspring's postnatal height growth seem to follow a differential pattern over time. There were stronger associations between the mother's height and the child's height growth velocity in the first months of life and stronger associations with the father's height after two years. If one admits that genes account for most of these associations, these mirroring patterns could reveal an underlying biological mechanism in growth regulation involving an alternating influence of maternal and paternal genes. Imprinting is well demonstrated in genes involved in the control of growth [START_REF] O'dell | Insulin-like growth factor II (IGF-II)[END_REF] and offers a potential mechanism behind our findings. For example, over-expression of paternal alleles in the IGF2 gene in Beckwith-Wiedemann syndrome which results in a marked acceleration of bone growth in infancy and early childhood [START_REF] Haig | Evolution in Health and Medicine Sackler Colloquium: Transfers and transitions: Parent-offspring conflict, genomic imprinting, and the evolution of human life history[END_REF] is in accordance with the stronger association between height growth velocity with father's than mother's height after two years in our results. A shared environment hypothesis (food intake, physical activity) between parents and offspring, enhancing the offspring's height growth during specific periods is less plausible.

The known relationship between offspring's ponderal index at birth and maternal but not paternal BMI is in accordance with the hypothesis of a preponderant impact of the maternal intra-uterine environment on fetal weight growth. The nutritional status of the mother before and during pregnancy, her metabolism and her dietary intake during pregnancy will determine the amount of glucose, oxygen and other nutrients the foetus will receive via the placenta [START_REF] Burk | Body size in five-year-old twins: heritability and comparison to singleton standards[END_REF][START_REF] Heude | Study of association between common variation in the insulin-like growth factor 2 gene and indices of obesity and body size in middle-aged men and women[END_REF].

In this context, genetic traits for weight growth, inherited from the mother or the father, may show up more frankly after birth [START_REF] Burk | Body size in five-year-old twins: heritability and comparison to singleton standards[END_REF][START_REF] Ong | Maternal-fetal interactions and birth order influence insulin variable number of tandem repeats allele class associations with head size at birth and childhood weight gain[END_REF]. In our study, paternal BMI was correlated with the offspring's early postnatal weight growth contrasting with the absence of association with birthweight. These relationships had already been highlighted in distinct studies [START_REF] Griffiths | Differential parental weight and height contributions to offspring birthweight and weight gain in infancy[END_REF][START_REF] Knight | Evidence of genetic regulation of fetal longitudinal growth[END_REF][START_REF] Vella | Association analysis of the IGF1 gene with childhood growth, IGF-1 concentrations and type 1 diabetes[END_REF] Consequently, paternal and maternal contributions to the offspring's BMI was similar at nine months. Unlike the study from Griffiths et al [START_REF] Griffiths | Differential parental weight and height contributions to offspring birthweight and weight gain in infancy[END_REF], the mother's BMI did not significantly contribute to the child weight growth velocity in the first three years of life. Our results concur with the observation that the relationship between body weight of mothers and their offspring's weight is weaker during the first two years of life than at birth [START_REF] Stunkard | Predictors of body size in the first 2 y of life: a high-risk study of human obesity[END_REF]. The postnatal preponderant influences of the fathers' BMI would unlikely be explained by environmental factors shared with the offspring during early infancy, but rather by paternal genes, which could not be expressed in utero for reasons such as maternal constraint or because they influence food intake behaviour. An interesting hypothesis has been advanced from the studies of human disorders of imprinted genes suggesting that alleles of paternal origin, expressed in infants, have been selected to favor more intense suckling than those of maternal origin [START_REF] Haig | Evolution in Health and Medicine Sackler Colloquium: Transfers and transitions: Parent-offspring conflict, genomic imprinting, and the evolution of human life history[END_REF]. In particular, children with Prader-Willi syndromecaused by defect of a paternally expressed alleleexperience a poor suck and failure to thrive in infancy. Therefore alleles of paternal origin may have a specific effect on offspring growth during the suckling periods which fits well with our observation of preferential association with paternal than maternal BMI in the fist months of life.

Maternal BMI was associated with the child's weight growth velocity from three years onwards. It is tempting to speculate that there is an increasing influence of shared environment and behaviours from three years onward. It has been shown that parents play a role in the eating behaviour of their offspring [START_REF] Fulkerson | Weight-related attitudes and behaviors of adolescent boys and girls who are encouraged to diet by their mothers[END_REF]. Preferences for a fatty diet have been related to parental adiposity [START_REF] Benton | Role of parents in the determination of the food preferences of children and the development of obesity[END_REF]. One could expect a preponderant role of the mothers, who usually spend more time with their children than fathers in shaping the eating behaviour of their child. However, a stronger association of the offspring's BMI or weight growth rate with the mother's BMI compared to the father's BMI is not consistently found [START_REF] Whitaker | The relationship of childhood adiposity to parent body mass index and eating behavior[END_REF]. Such discrepancies might be dependant on the relative involvement (eg time spent) of the father compared to the mother across the studies, and on the offspring's gender. In our study, the contribution of the mother's BMI was significantly higher in girls than in boys.

Our results about the association between parental BMI and the child's weight growth velocity can also be considered in the context of our previous finding of two periods where weight growth velocity predicts fat mass in adolescence [START_REF] Botton | Postnatal weight and height growth velocities at different ages between birth and 5 y and body composition in adolescent boys and girls[END_REF]. The first periodthe first six months of lifeis more related to paternal BMI, which suggests a role for genetic factors.

The second periodfrom three years onwardsis related to both parents' BMI. It is noteworthy that during the intermediate period where growth velocity was not related to fat mass in adolescence, there was also no significant contribution of the parents' BMI.

A first limitation of our study is that parental measurements were performed during the FLVS II study in 1999, when the offspring was already 8-17 years old. However, there is a high tracking of BMI during adulthood. Adjusting for the offspring's age at the date of parental measurements did not change the results. Obviously, this limitation does not concern relationships with height which, given the parents' age, should not change between their child's birth and the examination. Secondly, our sample, while being large enough to highlight relationships between parental and child characteristics, may not be powerful enough to highlight all differences between the mother-child and the father-child relationships, nor interactions with offspring's gender. However, stratification by gender allowed us to point out the significantly stronger contributions of the mother than the father to their daughter's weight and height growth velocities from five years onwards. Thirdly, we did not have data to test for a difference of the parental contributions according to whether the children were breastfed or not.

A first strength is the quality of the parental data, measured by trained clinicians, whereas other epidemiological studies often use reported data, particularly for the fathers. Secondly, the growth model allowed us to study the dynamic aspect of the relationships between parents' and children's characteristics: whereas measurements at given ages reflect the cumulative past growth experience, including the prenatal period, weight and height growth velocities are instantaneous estimates at given time points and enable identification of effects at specific periods.

In conclusion, we have shown successive periods from birth to seven years, during which the growth in weight and height have preferential and successive associations with either paternal or maternal body size. It is tempting to parallel this finding with the current knowledge on growth physiology in infancy, where growth hormone independent stimuli control growth in prenatal and early postnatal life, followed by a take-over by growth hormone dependent mechanisms [START_REF] Veldhuis | Endocrine control of body composition in infancy, childhood, and puberty[END_REF]. The recent hypotheses raised from the analysis of human disorders of imprinted genes offers another conceptual framework for the interpretation of our results [START_REF] Haig | Evolution in Health and Medicine Sackler Colloquium: Transfers and transitions: Parent-offspring conflict, genomic imprinting, and the evolution of human life history[END_REF].

Our results should trigger the search for specific genetic, epigenetic or shared environmental factors from the mother and the father that influence childhood growth and thus have major implications for further studies of height and weight growth regulation, including genetic studies.

FIGURE LEGENDS:

1A: Mean weight growth velocities (and SD) in boys and girls at different ages between birth and seven years.

1B:

  Mean height growth velocities (and SD) in boys and girls at different ages between birth and seven years. 2A: Coefficients of a multiple regression model for each age (and 95% confidence intervals) corresponding to the increase in offspring's ponderal index at birth or BMI at nine months and two years for a 1 kg/m 2 increase in father's or mother's BMI. Models included parent's BMIs, mother's age, offspring's gender, gestational age when model included birth measurements and took into account the family resemblance. 2B: Coefficients of a regression model for each selected age (and 95% confidence intervals) corresponding to the increase in offspring's weight growth velocity at different ages between 3 months and 7 years for a 1 kg/m 2 increase of in father's or mother's BMI. Models included parent's BMIs, mother's age, offspring's gender and took into account the family resemblance. 3A: Coefficients of a multiple regression model for each age (and 95% confidence intervals) corresponding to the increase in offspring's length at birth, nine months and 2 years for a 10 cm increase in father's or mother's height. Models included parent's heights, mother's age, offspring's gender, gestational age when model included birth parameters and took into account the family resemblance. 3B: Coefficients of a multiple regression model for each age (and 95% confidence intervals) corresponding to the increase in offspring's height growth velocity at different ages between 3 months and 7 years for a 10 cm increase in father's or mother's height. Models included parent's heights, mother's age, offspring's gender and took into account the family resemblance.

Acknowledgment:

We thank directors and teachers of the schools who made the study possible. We thank the parents and their children who agreed to participate in the study.

We also thanks Pr Yves Le Bouc, Dr S Cabrol, Dr I Netchine, INSERM unit 535 and Hôpital Armand Trousseau for their helpful comments on our results.

Members of the Fleurbaix Laventie Ville Santé Group: Arnaud Basdevant, Jean-Michel Borys, Jean-Louis Bresson, Marie-Aline Charles, Pierre Ducimetière, Philippe Froguel, Barbara Heude, Agnes Lommez, Jean-Michel Oppert, Monique Romon.

Marie-Aline Charles received grants from the "Association de Langue Française pour l'étude des diabètes et du metabolisme", ALFEDIAM; the "Mutuelle Generale de l'Education Nationale", MGEN; and the study was supported by Laboratoires Knoll, CEDUS, Groupe Fournier, Lesieur, Nestlé France, Produits Roche and TEPRAL. The study sponsors had no role in the study design, in the collection analysis and interpretation of data; in the writing of the manuscript, and in the decision to submit the manuscript for publication.