# Generating functions of stochastic L-systems and application to models of plant development.

If the interest of stochastic L-systems for plant growth simulation and visualization is broadly acknowledged, their full mathematical potential has not been taken advantage of. In this article, we show how to link stochastic L-systems to multitype branching processes, in order to characterize the probability distributions and moments of the numbers of organs in plant structure. Plant architectural development can be seen as the combination of two subprocesses driving the bud population dynamics, branching and differentiation. By writing the stochastic L-system associated to each subprocess, we get the generating function associated to the whole system by compounding the associated generating functions. The modelling of stochastic branching is classical, but to model differentiation, we introduce a new framework based on multivariate phase-type random vectors.Keywords: Stochastic L-system, multitype branching process, phase-type, generating functions, plant development, plant growth models, GreenLab

1 Introduction

Models of plant development (or organogenesis) describe the dynamic creation of organs (leaves, internodes, flowers/fruits) and how they arrange to form plant structures. When the smallest scale of interest is that of organs (and not cells), discrete models are generally used to simulate plant structural development. The parallel rewriting grammar introduced by Lindenmayer (1968) (called L-system) is particularly adapted to model the evolution of branching patterns and its algorithmic power has been broadly taken advantage of since Smith (1984). The stochastic version of this type of grammar gives interesting results from simulation and graphical points of view by increasing the realistic aspect of geometric plants, cf. Prusinkiewicz and Lindenmayer (1990). However, when the aim of development models is not purely geometric but a basis for the study of plant functioning (biomass production and repartition), there is a need to bypass the full simulation of plant structure and to write properly the dynamic equations of development, in order to compute the growth biophysical processes, cf. De Reffye et al. (2008). A solution based on grammar factorization was proposed by de Reffye et al. (2003) and Cournede et al. (2006). It was also derived in the stochastic case in Kang et al. (2007) and the formalism of stochastic grammars was taken advantage of to compute the generating functions of simple branching plant structures.

The objective of this article is to explore the relationship between stochastic L-Systems and multi-type branching processes in order to compute the corresponding probability generating functions and the distribution moments of the number of organs. The formalism is applied to a stochastic model of plant development (GreenLab). The article starts with the presentation of the basic botanical concepts and explains how a stochastic grammar can be derived. The theoretical framework is then introduced and we show how to derive generating functions from stochastic L-systems. It is applied to the two main botanical processes underlying plant development, branching and differentiation sequence of the meristem, which are first studied separately and then composed to form the complete stochastic model. The modelling of stochastic branching is classical but, to model differentiation, we introduce a new framework based on multivariate phase-type random vectors. Finally, the recursive relationship to compute the generating functions of plant structures can be introduced and applied to obtain the moments of the numbers of organs.

2 Botanical concepts to model plant development

As explained in Barthelemy and Caraglio (2007), organogenesis results from the functioning of undifferentiated cells constituting the apical meristem and located at the tip of axes. When in active phase, this meristem forms buds that will develop into new growth units composed of one or several metamers (also called phytomers). A metamer is a botanical entity chosen as the elementary scale to model plant architectural development in this study. It is composed of an internode bearing organs: axillary buds, leaves, flowers. Depending on species, metamers are set in place rhythmically--plants grow by successive shoots of several metamers--or continuously--when meristems keep on functioning and generate metamers one by one. In both cases, the time between the appearances of new shoots defines the architectural Growth Cycle. A Growth Unit is the set of metamers built by a bud during a growth cycle (only one in the continuous case). These metamers can be of different kinds and ordered according to botanical rules, like acrotony. For example, most temperate trees grow rhythmically, new shoots appearing at spring. If we do not consider polycyclism and neoformation, the architectural growth cycle corresponds to one year.

In this work, we do not consider time scales that are smaller than the architectural growth cycle and we study the development of new growth units as a discrete process. The Chronological Age (CA) of a plant (or of an organ) is defined as the number of growth cycles it has existed for.

Since metamers may bear axillary buds, plant architecture develops into a hierarchical branching system. Barthelemy and Caraglio (2007) underlined that architectural units can be grouped into categories characterized by a particular combination of morphological parameters. Thus, the concept of Physiological Age (PA) was introduced to represent the different types of growth units and axes. For instance, on coffee trees, there are two types: orthotropic trunk and plagiotropic branches. The main trunk's PA is equal to 1 and the oldest PA denoted by P corresponds to the ultimate state of differentiation for an axis, it is usually short, without branches. We need less than 5 PAs to describe the axis typology of most trees. The apical meristem or bud of an axis is thus characterized by the PA of the growth unit that it may produce and a metamer is characterized by its PA i (which is the PA of the growth unit that it belongs to). Moreover, along an axis, the morphological features of the growth unit may evolve with the age of the apical meristem. This process is described as the meristem sequence of differentiation by Barthelemy and Caraglio (2007), and corresponds to a transition to a superior PA of the meristem.

Since the number of potential axillary buds per metamer can be considered as a fixed botanical data, it is straightforward to deduce the whole plant structure from the population of buds. Moreover, we will often identify the meristem with the bud that it forms, for the sake of simplicity. Therefore, in the following, modelling plant development is equivalent to studying the dynamic evolution of the population of buds. It is mainly driven by the two botanical processes described above: branching resulting from the appearance of lateral buds in growth units and the differentiation sequence of meristems resulting in the change of PAs for terminal buds. As a consequence, at growth cycle n, a bud is characterized by 3 indices: its PA [phi], the CA k of the axis (which is also the CA of the meristem) and the initial PA of the meristem [beta]. It will be denoted by [b.sup.k,[beta].sub.n,[phi]]. In the sequel, the two botanical processes will be referred to as branching and differentiation.

[FIGURE 1 OMITTED]

The productions of buds deeply vary according to their positions, environmental conditions, plant age ... The deterministic modelling of such phenomena is difficult. For this reason, the development processes are considered as stochastic. For the sake of clarity, we only consider in this study growth units with one metamer. However, the same principles would apply for complex growth units with a random number of metamers as described in Kang et al. (2008). The computations of the general case are given in Loi and Cournede (2008). Experimental observations on plant populations and statistical analyses allow the estimation of the parameters of stochastic organogenesis models, cf. de Reffye (1979).

* [P.sub.l]: bud survival probability. At each growth cycle, a bud may stay alive with probability [P.sub.l] or die with probability (1 - [P.sub.l]). It may depend on bud's PA and will thus be denoted by [P.sub.l](i) for 1 [less than or equal to] i [less than or equal to] P.

* [P.sub.a]: bud activity probability. At each growth cycle, if a bud is alive, it may stay dormant with probability (1-[P.sub.a]) or produce a new growth unit with probability [P.sub.a]. [P.sub.a] may also depend on bud's PA: [P.sub.a](i) for 1 [less than or equal to] i [less than or equal to] P.

* [P.sup.b.sub.i,j](k): production probabilities. If at a given growth cycle, a bud of PA i is active, [P.sup.b.sub.i,j](k) is the probability that the growth unit it develops into bears k axillary buds of PA j. Botanical constraints usually impose that 0 [less than or equal to] k [less than or equal to] [B.sup.max.sub.i,j].

* [[lambda].sub.i]: inverses of average occupation times. During its sequence of differentiation, a meristem stays of PA i for an average

period of 1/[[lambda].sub.i].

* [q.sub.i,j] : transition probabilities. When a meristem of PA i changes its PA, [q.sub.i,j] is the probability that its new PA equals j. Note that the botanical differentiation sequence imposes that [q.sub.i,j] = 0 if j [less than or equal to] i.

In computational models, plants are generally represented as words in a language based on a generative parallel rewriting grammar also called L-system (Lindenmayer (1968), Smith (1984), Prusinkiewicz and Lindenmayer (1990), Francon (1990)). It is particularly adapted to describe the production of buds at discrete time steps. The following section will introduce a probabilistic framework to analyze stochastic models of plant development, based on stochastic L-systems.

3 A probabilistic framework based on stochastic L-systems

In this section, V = {[v.sub.1]; [v.sub.2], ..., [v.sub.m]} denotes an alphabet, W the set of all words over V and [W.sup.+] the set of nonempty words over V . Let 1 be the empty word and "." be the concatenation operator. Then, (W, ., 1) is a noncommutative monoid.

3.1 Markov kernel and stochastic L-system

We first recall some basic definitions and properties of the Markov chain theory (Stroock (2005)). Throughout the following all sets are finite or countable.

Definition 3.1 (Transition matrix) A transition matrix (or kernel) from A to B is a map k : A x B [right arrow] R such that k(a, b) [greater than or equal to] 0 and [[summation].sub.b[member of]B] k(a, b) = 1 for all a [member of] A. A transition matrix or Markov kernel on A is a transition matrix from A to A.

Definition 3.2 (Markov chain) A Markov chain with state space S is a stochastic process [([X.sub.n]).sub.n[member of]N] on some probability space ([OMEGA], F, P), such that P([X.sub.n+1] = y | [X.sub.n] = x) = [p.sub.x,y] for all x, y [member of] S. P = ([p.sub.x,y]) is a Markov kernel on S.

Theorem 3.1 For every tupel (P, [mu]) of Markov kernel P on S and probability measure [mu] on S there exists a Markov Chain [([X.sub.n]).sub.n[member of]N] of kernel P and starting measure [mu] = P([X.sub.0] [member of] x).

We propose beneath a definition for stochastic 0L-systems on V based on the notion of Markov kernel.

Definition 3.3 (Stochastic 0L-system) A stochastic 0L-system is a construct G = <[[omega].sub.a], [pi]> where:

* [[omega].sub.a] [member of] [W.sup.+] is called an axiom. It represents the structure initiating the growth.

* [pi] is a transition matrix from V to W.

Remark: Definition 3.3 is equivalent to the one given by Prusinkiewicz and Lindenmayer (1990): the set [P.sub.r] = {(x, y) [member of] V x W | [[pi].sub.x,y] > 0} defines the production rules. It represents the possible evolutions for all letters throughout the L-system.

We can now define a more general class of L-systems called stochastic F0L-system, extending the classical definition of F0L-system (Rozenberg and Salomaa (1980), p. 89) to the stochastic case:

Definition 3.4 (Stochastic F0L-system) A stochastic F0L-system is a construct G = <A, [pi]> where:

* A is a finite nonempty subset of V (called the set of axioms of G).

* for every [[omega].sub.a] [member of] A, G[[[omega].sub.a]] = <[[omega].sub.a], [pi]> is a stochastic 0L-system (called component system of G).

Proposition 3.1 Let P = ([p.sub.x,y]) be the square matrix on W such that, for all x = [x.sub.1][x.sub.2] ... [x.sub.n] [member of] W with ([x.sub.1], [x.sub.2], ..., [x.sub.n]) [member of] [V.sup.n] and for all y [member of] W:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Then, P is a Markov kernel on W. For a stochastic F0L-system <A, [pi]>, P is called the associated Markov kernel.

Remark: To every component system G[[[omega].sub.a]] = <[[omega].sub.a], [pi]> we associate a Markov chain [([F.sub.n] [[[omega].sub.a]]).sub.n[member of]N] via the associated kernel P of G and the starting measure [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] where [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] is the Dirac measure concentrated on [[omega].sub.a].

Beside plant topological structure directly deduced from L-Systems, in order to compute plant functioning, the numbers of organs are crucial variables (see for example De Reffye et al. (2008)). To determine them, the order of symbols in words does not play any role, and we can consider the L-systems as commutative. Let R be an equivalence relation on W defined as follows: [w.sub.1]R[w.sub.2] [??] there exists [PI], a permutation on the symbol ranks, such that [PI]([w.sub.1]) = [w.sub.2]. Let us denote the quotient set W/R by [W.sup.*]. From now on, each word w [member of] W will be assimilated to the ordered representative [w.sup.*] of its equivalence class (i.e. wR[w.sup.*] and [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] for ([[alpha].sub.1], ..., [[alpha].sub.m]) [member of] [N.sup.m]). [W.sup.*] is isomorphic to [N.sup.m]. Let Y be the canonical isomorphism from [W.sup.*] into [N.sup.m]. For a word [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] in [W.sup.*], we have Y(w) = ([[alpha].sub.1], ..., [[alpha].sub.m]). Let us denote by [Y.sub.i](w) the i-th component of Y(w) (i.e. [Y.sub.i](w) = [[alpha].sub.i]).

In the following, the transition matrix [pi] is thus considered as a map from V x [W.sup.*] into R.

3.2 Multitype branching processes and stochastic L-systems

Stochastic L-systems are closely related to multitype branching processes. Let us recall first the definition of a multitype branching process (Harris (1963), Athreya and Ney (2004)). Let us consider a population with m types of individuals. Assume a type i individual produces children of all types according to a probability distribution {[P.sub.i](j) : j = ([j.sub.1], ..., [j.sub.m]), [j.sub.i] [member of] N, 1 [less than or equal to] i [less than or equal to] m}. Assume all individuals produce offspring independently of each other and of the past history of the process. Let [Z.sub.n,i] be the number of type i individuals in the n-th generation. Let {[[xi].sup.(k).sub.n,i] : n [member of] [N.sup.*], k [member of] N, 1 [less than or equal to] i [less than or equal to] m} be independent random vectors in [N.sup.m] with [[xi].sup.(k).sub.n,i] having distribution [P.sub.i](.).

Definition 3.5 (multitype Galton-Watson branching process) If the vector [Z.sub.n] = ([Z.sub.n,1], ..., [Z.sub.n,m]) of population sizes in the n-th generation evolves by the recursive relation

[Z.sub.n+1] = [m.summation over (i=1)] [[Z.sub.n.i].summation over (k=1)] [[xi].sup.(k).sub.n,i], (1)

then [([Z.sub.n]).sub.n[member of]N] is a multitype Galton-Watson branching process.

The j-th component of [[xi].sup.(k).sub.n,i] represents the number of type j individuals produced by the k-th type i individual in the n-th generation. The set [{[P.sub.i](.)}.sub.i[member of]{1, ..., m}] is called the offspring distribution. Considering this definition, we have the following theorem:

Theorem 3.2 Let G = <A, [pi]> be a stochastic F0L-system on V = {[v.sub.1], ..., [v.sub.m]} with letters acting independently. Let G[[[omega].sub.a]] = <[[omega].sub.a], [pi]> be a component system of G and [([F.sub.n] [[[omega].sub.a]]).sub.n[member of]N] the corresponding Markov chain. Then, [(Y([F.sub.n] [[[omega].sub.a]])).sub.n[member of]N] is a multitype Galton-Watson branching process whose offspring distributions [{[P.sub.i](.)}.sub.i[member of]{1, ..., m}] are given by:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Proof: Let [{[P.sub.i](.)}.sub.i[member of]{1, ..., m}] be a set of probability distributions defined as follows:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Let {[[xi].sup.(k).sub.n,i] : n [member of] [N.sup.*] k [member of] N, 1 [less than or equal to] i [less than or equal to] m} be independent random vectors with [[xi].sup.(k).sub.n,i] having distribution [P.sub.i](.). Thus, the q-th component of [[xi].sup.(k).sub.n,i] is the random variable that gives the number of type "[v.sub.q]" symbols in a word generated by a type "[v.sub.i]" symbol throughout G. Since there are [Y.sub.i] ([F.sub.n][[[omega].sub.a]]) type "[v.sub.i]" symbols in [F.sub.n][[[omega].sub.a]], we have:

[Y.sub.q] ([F.sub.n+1][[[omega].sub.a]]) = [m.summation over (i=1)] [[Y.sub.i] ([F.sub.n][[[omega].sub.a]]).summation over (k=1)] [([[xi].sup.(k).sub.n,i]).sub.q]

and then,

Y([F.sub.n+1][[[omega].sub.a]]) = [m.summation over (i=1)] [[Y.sub.i]([F.sub.n][[[omega].sub.a]]).summation over (k=1)] [[xi].sup.(k).sub.n,i]

N.B.: the converse is true. Giving a multitype Galton-Watson branching process, there is a stochastic L-system whose stochastic production function is entirely determined by the offspring distribution.

Let us now define the generating function associated to a stochastic 0L-system. Let S = ([s.sub.1], ..., [s.sub.m]) [member of] [[0, 1].sup.m].

Definition 3.6 (generating function associated to a stochastic 0L-system) Let G = <A, [pi]> be a stochastic F0L-system on V = {[v.sub.1], ..., [v.sub.m]} with letters acting independently. Let G[[[omega].sub.a]] = <[[omega].sub.a], [pi]> be a component system of G and [([F.sub.n][[[omega].sub.a]]).sub.n[member of]N] the corresponding Markov chain. For n [member of] N, the generating function [[psi].sub.n][[[omega].sub.a]] of [F.sub.n][[[omega].sub.a]] is defined as follows:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

[[psi].sub.1] [[[omega].sub.a]] is said to be the generating function associated to G[[[omega].sub.a]].

By using the classical composition of generating functions (Harris (1963)) for a multitype Galton-Watson branching process, we deduce directly the following theorem:

Theorem 3.3 Let G = <V, [pi]> be a stochastic F0L-system on V = {[v.sub.1], ..., [v.sub.m]} with letters acting independently. For all v [member of] V, let G[v] = <v, [pi]> be a component system of G and and [([F.sub.n][v]).sub.n[member of]N] the corresponding Markov chain. For all n [member of] N, let [[psi].sub.n][v] be the generating function associated to [F.sub.n][v].

Then,

[for all]n [member of] N, [for all][[omega].sub.a] [member of] V, [[psi].sub.n+1][[[omega].sub.a]](S) = [[psi].sub.1] [[[omega].sub.a]] ([[psi].sub.n] [[v.sub.1]](S), ..., [[psi].sub.n][[v.sub.m]](S))

4 Application to stochastic plant development

The GreenLab model of plant development is based on the botanical modelling concepts recalled in section 2. It can be described by a stochastic F0L-system [G.sup.tot] = <V, [[pi].sup.tot]> (Kang et al. (2007)). Let N [member of] N be the time during which we observe the growth of the plant. A bud [b.sup.k,[beta].sub.n,[phi]] is symbolized by [s.sup.k,[beta].sub.[phi]] and a metamer of PA j by [m.sub.j]. V is the union of S = {[s.sup.k,[beta].sub.[phi]] : ([phi], [beta]) [member of] [{0, ..., P}.sup.2], k [member of] {1, ..., N}} and M={[m.sub.j] : j [member of] {1, ..., P}} where S is the set of buds (nonterminal elements) and M is the set of metamers (terminal elements). A dead bud will be represented by the empty word 1.

The seed can be considered as the initial bud and is thus represented by [s.sup.0,1.sub.1]. We consider the component system [G.sup.tot][[s.sup.0,1.sub.1]] and [([F.sup.tot.sub.n] [[s.sup.0,1.sub.1]]).sub.n[member of]N] the corresponding Markov chain. For n [member of] N, [F.sup.tot.sub.n] [[s.sup.0,1.sub.1]] is the random variable representing the possible realisations of a plant after n growth cycles. Our objective is to determine the generating function associated to [F.sup.tot.sub.n] [[s.sup.0,1.sub.1]] with Theorem 3.3. Thus, we only need to determine the generating function associated to [F.sup.tot.sub.1] [s.sup.0,1.sub.1]]. However, the development processes underlying GreenLab organogenesis are complex. The dynamics of the population of buds results from the combination of both branching and differentiation (see Section 2) and we do not have an easy access to the distribution of [F.sup.tot.sub.1] [[s.sup.0,1.sub.1]]. Thus, the idea is to break down the whole system (branching + differentiation) into two subsystems and study them separately. This is the aim of Sections 4.1 (for branching) and 4.2 (for differentiation). In both cases, we write the corresponding generating function. Section 4.3 combines the two processes by using the results thus obtained and gives the recursive equations for the generating functions derivating from [G.sup.tot]. The recursive equations for the moments of the numbers of metamers are finally deduced.

4.1 Branching

In this section we focus on the GreenLab development model with only branching. The differentiation process is not taken into account. This section gives not only the generating function for the branching process but also a method to build it. The model can be represented by a stochastic F0L-system [G.sup.br] = <V, [[pi].sup.br]>. V is defined as in Section 4. For every component system [G.sup.br][v] of [G.sup.br], the corresponding Markov chain is denoted [(F.sup.br.sub.n] [v]).sub.n[member of]N. In order to illustrate the branching process, let us take the example of coffee trees, plants with stochastic branching and no differentiation. They have two PAs (i.e. P=2). A metamer of PA 1 can bear a maximum of two lateral buds of PA 2 (i.e. [B.sup.max.sub.1,2] = 2). A metamer of PA 2 does not bear any lateral bud. At each growth cycle, there are three possible evolutions for a bud:

* case 1: [F.sup.br.sub.1] [[s.sup.k,[phi].sub.[phi]] = 1 with [phi] [member of] {1, 2}, the bud dies.

* case 2: [F.sup.br.sub.1] [[s.sup.k,[phi].sub.[phi]] = [s.sup.k+1,[phi].sub.[phi]] with [phi] [member of] {1, 2}, the bud is still alive but rests, .

* case 3: the bud is active. In that case, it produces a new metamer with lateral buds ([F.sup.br.sub.1] [[s.sup.k,1.sub.1]] = [m.sub.1][s.sup.k+1,1.sub.1] [([s.sup.0,2.sub.2]).sup.j] with j [member of] {1, 2}) or not ([F.sup.br.sub.1] [[s.sup.k,[phi].sub.[phi]]] = [m.sub.[phi]][s.sup.k+1,[phi].sub.[phi]] with [phi] [member of] {1, 2}).

A probability of occurrence can be given for each of these possible evolutions. Let us now consider the general case. The aim is to write the generating function associated to [G.sup.br]. Let [[psi].sup.br.sub.1] [[s.sup.k,[beta].sub.[phi]]](S,M) be the generating function associated to the component system [G.sup.br][[s.sup.k,[beta].sub.[phi]]. S = [([s.sup.k',[beta]'.sub.[phi]]).sub.k', [beta]', [phi]'] and M = [(m.sub.[phi]').sub.[phi]'] are vectors respectively on [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] and [[0, 1].sup.P]. Then,

[FIGURE 2 OMITTED]

Theorem 4.1 The generating function associated to [G.sup.br][[s.sup.k,[beta].sub.[phi]]] is given by

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Proof: Using Definition 3.6, we get the generating function by enumerating all the possible evolutions of [s.sup.k,[beta].sub.[phi]].

* case 1: [F.sup.br.sub.1] [[s.sup.k,[beta].sub.[phi]]] = 1 and the probability of occurrence is 1 - [P.sub.l]([phi]).

* case 2: [F.sup.br.sub.1] [[s.sup.k,[beta].sub.[phi]]] = [s.sup.k+1,[beta].sub.[phi]] and the probability of occurrence is [P.sub.l]([phi])(1 - [P.sub.a]([phi])).

* case 3: the rules are of type [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] with probability of occurrence

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Therefore, by summing the three cases, we get the equation of Theorem 4.1. 2

4.2 Differentiation

In this section we focus on the differentiation process, corresponding to the change of meristem's PA. The branching process is not taken into account. We proceed as in Section 4.1 and write the stochastic F0L-System [G.sup.dif] = <V, [[pi].sup.dif]> associated to the differentiation process in order to write the corresponding generating function. V is defined as in Section 4. For every component system [G.sup.dif] [v] of [G.sup.dif], the corresponding Markov chain is denoted [(F.sup.dif.sub.n] [v]).sub.n[member of]N].

To illustrate differentiation, let us consider the example of maize. In standard cultivation conditions, maize is a mono-stem plant, that is to say without ramification. However, we can distinguish two types of metamers along the stem. The first ones are short and can potentially bear tillers. They are followed by longer phytomers after meristem differentiation. Finally, the meristem ends up by flowering, which terminates the differentiation sequence and the stem development, see Figure 3. The two types of metamers are characterized by two different PAs (P = 2). We assume that the terminal bud of an axis of CA k = 0 can not change its PA immediately. At each growth cycle, there are three possible evolutions for a bud through the differentiation process:

* case 1: [F.sup.dif.sub.1] [[s.sup.k,1.sub.[phi]]] = [s.sup.k,1.sub.[phi]] with [phi] [member of] {1, 2}, the bud does not differentiate. The PA of the terminal bud remains the same.

* case 2: [F.sup.dif.sub.1] [[s.sup.k,1.sub.1]] = [s.sup.k,1.sub.2], the bud differentiates and is still alive. The PA of the terminal bud changes and is higher (cf Section 2)

* case 3: [F.sup.dif.sub.1] [[s.sup.k,1.sub.[phi]] = 1 with [phi] [member of] {1, 2}, the bud differentiates and dies. The growth of the axis stops.

[FIGURE 3 OMITTED]

A probability of occurrence can be given for each of these evolutions. Let us now consider the general case. In order to write the generating function, we first have to identify the stochastic process underlying the differentiation of the terminal bud.

The aim is to study the evolution of the PA [phi] of a bud with respect to the CA k of the axis that bears the bud. Let [b.sup.k,[beta].sub.n,[phi]] be a bud of PA [phi] at growth cycle n. The axis that bears the bud is of CA k and was initiated by a bud of PA [beta]. Then, we assume [phi] = [phi](k). Moreover, we assume that the occupation time of a bud in the PA j follows an exponential law of parameter [[lambda].sub.j] where [[lambda].sub.j] is defined in Section 2. It is in keeping with the botanical theory of meristematic differentiation sequence as a sequence of events inducing plant morphological ageing. Therefore, we can represent the discrete process [([phi](k)).sub.k[member of]N]] by a finite state space Markov process [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII].

We choose to model the stochastic behaviour of [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] with multivariate phase-type random vectors, which are particularly used in reliability theory (Neuts (1975)). An interesting aspect of phase-type distributions is that they can be written in a closed form, and consequently, various quantities of interest can be evaluated with relative ease (Assaf et al. (1984)). Let [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] be a right-continuous Markov process on a finite state space E = {1, ..., P, [DELTA]}. Assume 1, ... P are transient and [DELTA] is absorbing.

Then, the infinitesimal generator Q of size P + 1 by P + 1 has the form

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

where A is a P x P matrix and e is a vector of size P with all its components set to 1.

Let [[GAMMA].sub.1], ..., [[GAMMA].sub.n] be nonempty stochastically closed subsets of E for X. A subset [GAMMA] of E is said to be stochastically closed for X if once X enters [GAMMA], it never leaves it. Let a be the starting measure on E with a([DELTA]) = 0. Let [alpha] be the vector of size P defined by [alpha] = (a(1), ..., a(P)). We define [T.sub.k] = inf{j, X(j) [member of] [[GAMMA].sub.k]}.

Definition 4.1 (multivariate phase-type random vector) If the following hypotheses are true:

1--[[intersection].sup.n.sub.k=1] [[GAMMA].sub.k] = [DELTA],

2--the absorption into [DELTA] is certain,

3--P([T.sub.1] > 0, ..., [T.sub.n] > 0) = 1, then ([T.sub.1], ..., [T.sub.n]) is called a multivariate phase-type random vector with representation ([alpha], A).

Considering now the Markov process [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] we have E = {1, ..., P, d} where 1, ..., P represent the possible PAs (i.e. transient states) and d represents a dead bud, or the terminal flowering state as in Maize, (i.e. an absorbing state). We assume [phi](0) [not equal to] d. In the sequel, it will be important to specify the initial state of the Markov process. Let [([[phi].sup.[beta]](t)).sub.t[member of]R] be the PA process such that [[phi].sup.[beta]](0) = [beta] with [beta] [member of] {1, ..., P}. Let [a.sup.[beta]] be the starting measure on E for [[phi].sup.[beta]]. We assume that [a.sup.[beta]](d) = 0. Let [[alpha].sup.[beta]] be the vector of size P defined by [[alpha].sup.[beta]] = ([[alpha].sup.[beta]](1), ..., [[alpha].sup.[beta]](P)). Given that [[phi].sup.[beta]](0) = [beta], we have [[alpha].sub.[beta]](j) = [[delta].sub.[beta]](j) where [[delta].sub.[beta]] is the Dirac measure concentrated on [beta]. The infinitesimal generator Q is given by [{[[lambda].sub.j]}.sub.j[member of] {1, ..., P}] and the transition probabilities of the underlying Markov chain [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. Let us recall that [q.sub.i,j] = 0 for j [less than or equal to] i (Section 2). Then,

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Let [{[[GAMMA].sub.k]}.sub.k[member of]{1, ..., P+1}] be a family of subsets of E defined as follows:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Given [q.sub.i,j] = 0 for j [less than or equal to] i, [{[[GAMMA].sub.k]}.sub.k[member of]{1, ..., P+1}] are stochastically closed subsets of E for [[phi].sup.[beta]]. [T.sup.[beta].sub.k] denotes the time when [[phi].sup.[beta]] enters [[GAMMA].sub.k]:

[T.sup.[beta].sub.k] = inf{j, [[PHI].sup.[beta]](j) [member of] [[GAMMA].sub.k]}]

Then, we have the following proposition:

Proposition 4.1 For all [beta] [member of] {1, ..., P}, ([T.sup.[beta].sub.[beta]+1]], ..., [T.sup.[beta].sub.P+1]]) is a multivariate phase-type random vector with representation ([[alpha].sup.[beta]], A]). Moreover, we have 0 < [T.sup.[beta].sub.[beta]+1]] < [T.sup.[beta].sub.[beta]+2]] < ... < [T.sup.[beta].sub.P+1]].

Proof: For all [beta] [member of] {1, ..., P}, we have to prove the three hypotheses:

1--it is obvious that [[intersection].sup.P+1.sub.k=[beta]+1] [[GAMMA].sub.k] = d.

2--according to Neuts (1975), the absorption into d is certain if and only if A is non singular. This is true because det [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

3--given the only nonzero component of [[alpha].sup.[beta]] is the [beta]-th one, [T.sup.[beta].sub.k] > 0 is always true for k [member of] {[beta] + 1, ..., P + 1}. Then, the third hypothesis is true.

According to Definition 4.1, ([T.sup.[beta].sub.[beta]+1]], ..., [T.sup.[beta].sub.P+1]]) is a multivariate phase-type random vector with representation ([[alpha].sup.[beta]], A). The second part of the proposition is obvious because [[GAMMA].sub.P+1] [subset] [[GAMMA].sub.P] [subset] ... [subset] [[GAMMA].sub.1].

We will now determine the generating function derivating from [G.sup.dif]. The probability distribution of [F.sup.dif.sub.1] [[s.sup.k,[beta].sub.[phi]]] is given naturally by the multivariate phase-type distributions.

Let [[psi].sup.dif.sub.1] [[s.sup.k,[beta].sub.[phi]](S) be the generating function associated to [G.sup.dif] [[s.sup.k,[beta].sub.[phi]]. S = [([s.sup.k', [beta]'.sub.[phi]']).sub.k',[beta],[phi]'] is a vector on [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. Let D be the event {[T.sup.[beta].sub.[phi]] < k - 1, [T.sup.[beta].sub.[phi]+1] > k - 1}. Let [g.sub.k] be a diagonal matrix of size P by P whose i-th diagonal element is equal to 1 if i [member of] [[GAMMA].sup.c.sub.k] and 0 otherwise. [[GAMMA].sup.c.sub.k] denotes the complement of [[GAMMA].sub.k] in E. Then,

Theorem 4.2 The generating function associated to [G.sup.dif] [[s.sup.k,[beta].sub.[phi]]] is given by

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

if 1 [less than or equal to] k [less than or equal to] N - 1 and by

[[psi].sup.dif.sub.1] [[s.sup.0,[beta].sub.[beta]]](S) = [s.sup.0,[beta].sub.[beta]]

if k = 0 and all the probabilities are explicit functions of the matrix A:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Proof: We will proceed as in Section 4.1. For the bud [s.sup.k,[beta].sub.[phi]], we know that [phi] enters [[GAMMA].sub.[phi]] before k - 1 but not [[GAMMA].sub.[phi]+1]. It explains why all the probabilities in Theorem 4.2 have the form P(x|D). There are 3 cases for bud differentiation:

* [F.sup.dif.sub.1] [[s.sup.k,[beta].sub.[phi]]] = [s.sup.k,[beta].sub.[phi]]: the bud does not change its PA. This is the case if [[phi].sup.[beta]] does not enter [[GAMMA].sub.k+1] before k. Therefore, the probability of occurrence is P([T.sup.[beta].sub.[phi]+1] > k + 1|D).

* [F.sup.dif.sub.1] [[s.sup.k,[beta].sub.[phi]]] = [s.sup.k,[beta].sub.[phi]'] with [phi]' > [phi]: the bud changes its physiological for [phi]'. It means [[phi].sup.[beta]] enters [[GAMMA].sub.[phi]'] before k but not [[GAMMA].sub.[phi]'+1]. Therefore, the probability of occurrence is P([T.sup.[beta].sub.[phi]'] < k + 1, [T.sup.[beta].sub.[phi]'+1] > k + 1|D).

* [F.sup.dif.sub.1] [[s.sup.k,[beta].sub.[phi]]] = 1 : the bud dies. It means that [[phi].sup.[beta]] enters [[GAMMA].sub.P+1] before k. Therefore, the probability of occurrence is P([T.sup.[beta].sub.P+1] < k|D).

By summing the three cases, we have the expression of [[psi].sup.dif.sub.1] [[s.sup.k,[beta].sub.[phi]]](S). The case of [s.sup.0,[beta].sub.[beta]] is particular because we assume that a bud that has just initiated a new axis can not change its PA.

Modelling differentiation with multivariate phase-type random vectors is interesting since all probabilities are given explicitly. It was proved in Assaf et al. (1984) that for all [t.sub.P+1] [greater than or equal to] [t.sub.P] [greater than or equal to] ... [greater than or equal to] [t.sub.[beta]+1] [greater than or equal to] 0,

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Thus, given [T.sup.[beta].sub.[beta]+1] < [T.sup.[beta].sub.[beta]+2] < ... < [T.sup.[beta].sub.P+1], we get the expression of the probabilities of Theorem 4.2.

4.3 The complete model

In this section, we establish the recursion formulas for the generating functions derivating from the complete L-system [G.sup.tot] = <V, [[pi].sup.tot]> defined at the beginning of Section 4. The aim is to write the recursion formulas considering [G.sup.tot] as a composition of the L-systems [G.sup.br] and [G.sup.dif]. Let us define first a compound L-system:

Definition 4.2 (compound L-system) Let G = <V, [pi]>, [G.sup.1] = <V, [[pi].sup.1]> and [G.sup.2] = <V, [[pi].sup.2]> be stochastic L-systems on the same alphabet V. Let P, [P.sup.1] and [P.sup.2] be the Markov kernel associated respectively to G, [G.sup.1] and [G.sup.2]. G is said to be the composition of [G.sup.1] by [G.sup.2] if P = [P.sup.1].[P.sup.2]. We write G = [G.sup.1] o [G.sup.2].

The complete GreenLab development model mixes branching and differentiation. We assume the priority is given to the branching process at each growth cycle. Then, each step of [G.sup.tot] begins with a step of [G.sup.br] and goes on with a step of [G.sup.dif]. It is thus obvious that [G.sup.tot] = [G.sup.br] o [G.sup.dif]. This result leads to the following theorem:

Theorem 4.3 For all MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] and M = [([m.sub.[phi]']).sub.[phi]'] [member of] [[0,1].sup.P],

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Proof: Since [G.sup.tot] = [G.sup.br] o [G.sup.dif], the generating function derivating from [G.sup.tot] is the compound function of the generating function derivating from [G.sup.dif] by that of [G.sup.br] i.e.:

[[psi].sup.tot.sub.1] [s.sup.k,[beta].sub.[phi]]] (S,M) = [[psi].sup.br.sub.1] [[s.sup.k,[beta].sub.[phi]]] ([[psi].sup.dif](S), M)

Then, using Theorem 3.3, we get the result.

Setting [[psi].sup.tot.sub.1] [[s.sup.k,[beta].sub.[phi]] (S,M) = 0 for k + n > N + 1 and [[psi].sup.br] = [([[psi].sup.br.sub.1] [[s.sup.k,[beta].sub.[phi]]]).sub.k,[beta],[phi]], we have:

Theorem 4.4 For all n [less than or equal to] N - 1, [[psi].sup.tot.sub.n+1](S,M) = [[psi].sup.br]([[psi].sup.dif]([[psi].sup.tot.sub.n](S,M)),M).

Let [m.sub.j] [[s.sup.k,[beta].sub.[phi]], n] be the number of type j metamers in a structure initiated by [s.sup.k,[beta].sub.[phi]] after n growth cycles and let [M.sub.n] be the matrix of size (N + 1)[P.sup.2] by P whose j-th column vector is [(E[[m.sub.j] [[s.sup.k,[beta].sub.[phi]],n]]).sub.k,[beta],[phi]]. Let [e.sub.k] be the vector of size k with all its components set to 1. We deduce the fundamental recursion equation for the expectations of the numbers of organs on each type of structures. It is important to note that the first line of the matrix corresponds to those of the whole plant.

Theorem 4.5 for N > 0,

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Proof: Using the properties of generating functions, we have [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. Derivating the equation from Theorem 4.4 and taking [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] we get Theorem 4.5. 2

In the same way, we can give recursion formulas for the variance. The detailed results, as well as examples of simulations and biological discussions, are given in Loi and Cournede (2008).

5 Conclusion

While the interest of stochastic L-systems for plant growth simulation and visualization is broadly acknowledged, its full mathematical potential to characterize the probability distributions and moments of the numbers of organs in plant structure has not been taken advantage of. The need for an analytic computation of these distributions is crucial in stochastic functional-structural models in order to derive the distribution or the moments of biomass production (cf. Kang et al. (2008) for preliminary results). It has led us to clearly formalize the link between stochastic L-systems and multi-type branching processes, and thus to derive an inductive relationship to compute the associated generating functions. This framework was applied successfully to the GreenLab organogenesis model, by decomposing the development process into two botanical sub-processes, branching and meristem differentiation. For the latter, multivariate phase type random vectors were introduced to describe the stochastic sequence of meristem differentiation. From the inductive relationship on the generating functions of the numbers of organs, we can determine the inductive relationship giving the moments of the distributions.

The next step is the link of stochastic organogenesis and functioning models. It is possible to derive the approximate moments of biomass production thanks to differential statistics. An interesting issue concerns the modelling of the probability distributions governing organogenesis as functions of biomass production or groth rate (cf. Mathieu et al. (2006)). The stochastic processes resulting from these interactions between development and functioning should lead us to improve the proposed formalism.

References

D. Assaf, N. Langberg, T. Savits, and M. Shaked. Multivariate phase-type distributions. Operations Research, 32(3):688-702, 1984.

K. Athreya and P. Ney. Branching Processes. Dover Publications, 2004.

D. Barthelemy and Y. Caraglio. Plant architecture: a dynamic, multilevel and comprehensive approach to plant form, structure and ontogeny. Annals of Botany, 99(3):375-407, 2007.

P.-H. Cournede, M.-Z. Kang, A. Mathieu, J.-F. Barczi, H.-P. Yan, B.-G. Hu, and P. de Reffye. Structural Factorization of Plants to Compute their Functional and Architectural Growth. Simulation, 82(7):427-438, 2006.

P. de Reffye. Modelisation de l'architecture des arbres par des processus stochastiques. Simulation spatiale des modeles tropicaux sous l'effet de la pesanteur. Application au Coffea robusta. PhD thesis, Universite Paris-Sud, Centre d'Orsay, 1979.

P. de Reffye, M. Goursat, J. Quadrat, and B. Hu. The Dynamic Equations of the Tree Morphogenesis Greenlab Model. Technical Report 4877, INRIA, 2003.

P. De Reffye, E. Heuvelink, D. Barthelemy, and P.-H. Cournede. Biological and mathematical concepts for modelling plant growth and architecture. In Encyclopedia of Ecology. Jorgensen, D.E. and Fath, B., elsevier edition, 2008.

J. Francon. Sur la modelisation informatique de l'architecture et du developpement des vegetaux. In 2eme Colloque International: L'Arbre. Institut de Botanique, Montpellier, France, 1990.

T. Harris. The theory of branching processes. Springer, Berlin, 1963.

M.-Z. Kang, P.-H. Cournede, J.-P. Quadrat, and P. de Reffye. A stochastic language for plant topology. In T. Fourcaud and X. Zhang, editors, Plant growth Modeling, simulation, visualization and their Applications. IEEE Computer Society (Los Alamitos, California), 2007.

M.-Z. Kang, P.-H. Cournede, P. de Reffye, D. Auclair, and B.-G. Hu. Analytical study of a stochastic plant growth model: application to the greenlab model. Mathematics and Computers in Simulation, 78 (1):57-75, 2008.

A. Lindenmayer. Mathematical models for cellular interactions in development. i. filaments with one-sided inputs. Journal of Theoretical Biology, 18:280-289, 1968.

C. Loi and P.-H. Cournede. Description of the GreenLab development model with stochastic L-systems and Monte-Carlo simulations. Technical report, INRIA, 2008.

A. Mathieu, P.-H. Cournede, and P. de Reffye. A dynamical model of plant growth with full retroaction between organogenesis and photosynthesis. ARIMA Journal, 4:101-107, 2006. URL http://www-direction.inria.fr/international/arima/004/00407.htm.

M. F. Neuts. Probability distributions of phase type. Liber Amicorum Prof. Emeritus H. Florin, pages 176-206. University of Louvain, Belgium, 1975.

P. Prusinkiewicz and A. Lindenmayer. The Algorithmic Beauty of Plants. Springer-Verlag, New-York, 1990.

G. Rozenberg and A. Salomaa. The Mathematical Theory of L-systems. Academic Press, New York, 1980.

A. Smith. Plants, fractals and formal languages. Computer Graphics (SIGGRAPH 84 Conference Proceedings), 18(3):1-10, 1984.

D. Stroock. An introduction to Markov processes. Springer, 2005.

Cedric Loi (1,2), Paul Henry Cournede (1,2)

(1) Ecole Centrale Paris, Laboratory of Applied Mathematics and Systems, 92290 Chatenay Malabry, France

(2) INRIA Saclay-Ile-de-France, DIGIPLANTE, 91400 Orsay, France

Printer friendly Cite/link Email Feedback | |

Author: | Loi, Cedric; Cournede, Paul Henry |
---|---|

Publication: | DMTCS Proceedings |

Article Type: | Report |

Geographic Code: | 4EUFR |

Date: | Jan 1, 2008 |

Words: | 7570 |

Previous Article: | Convergence to the coalescent and its relation to the time back to the most recent common ancestor. |

Next Article: | A functional limit law for the profile of plane-oriented recursive trees. |

Topics: |