Abstract
We review and extend the formalism introduced by Peliti, that maps a Markov process to a path-integral representation. After developing the mapping, we apply it to some illustrative examples: the simple decay process, the birth-and-death process, and the Malthus-Verhulst process. In the first two cases we show how to obtain the exact probability generating function using the path integral. We show how to implement a diagrammatic perturbation theory for processes that do not admit an exact solution. Analysis of a set of coupled Malthus-Verhulst processes on a lattice leads, in the continuum limit, to a field theory for directed percolation and allied models.
Path integrals and perturbation theory for stochastic processes
Ronald Dickman; Ronaldo Vidigal
Departamento de Física, ICEx, Universidade Federal de Minas Gerais, 30123-970, Belo Horizonte, MG, Brasil
Address to correspondence Address to correspondence Ronald Dickman, Ronaldo Vidigal Departamento de Física, ICEx, Universidade Federal de Minas Gerais, 30123-970, Belo Horizonte, MG, Brasil
ABSTRACT
We review and extend the formalism introduced by Peliti, that maps a Markov process to a path-integral representation. After developing the mapping, we apply it to some illustrative examples: the simple decay process, the birth-and-death process, and the Malthus-Verhulst process. In the first two cases we show how to obtain the exact probability generating function using the path integral. We show how to implement a diagrammatic perturbation theory for processes that do not admit an exact solution. Analysis of a set of coupled Malthus-Verhulst processes on a lattice leads, in the continuum limit, to a field theory for directed percolation and allied models.
1 Introduction
It is often noted that nonequilibrium statistical mechanics lacks the comprehensive formalism of ensembles that has proved so useful in equilibrium. The reason is that equilibrium statistical mechanics treats stationary states for a special class of systems, that possess detailed balance. This allows one to bypass the dynamics, and study stationary properties directly. For systems out of equilibrium, we generally do not have such a shortcut, and must deal with the full dynamical problem, even if our goal is only to obtain stationary properties. In the case of stochastic systems with a discrete state space, the fundamental description is given by the master equation, which governs the evolution of the probability distribution. This class of problems includes a wide range of systems of current interest, that exhibit phase transitions or scale invariance far from equilibrium: driven lattice gases, birth-and-death processes such as directed percolation or the contact process, sandpile models, and interface growth models.
One of the more powerful tools for studying stochastic models is a formalism that maps the process to a path-integral representation. This mapping generates an effective action that can be studied using the tools of equilibrium statistical physics, for example, the renormalization group. Several methods for mapping a stochastic process to an equilibrium-like action have been proposed [1-7]. In this article we review the method developed by Peliti [8], and apply it to some simple stochastic processes. This method has several advantages. With it, one can map any birth-and-death process to a path-integral representation without ambiguity. In particular, the step of writing a Langevin equation, and of postulating noise autocorrelations, does not arise in this formalism. Thus it provides a direct path from the model of interest to an effective action, and (in the continuum limit), to its field theory, without the uncertainties that often attend the specification of the noise term [9]. A second advantage, which we explore in detail, is that it leads to a systematic perturbative analysis for Markov processes.
The principal aim of this article is to acquaint the reader with the formalism and provide a set of worked examples whose mastery will allow one to apply the method to problems at the frontier of research. While Peliti's article [8] provides an excellent exposition of the mapping, we include, for completeness, a derivation of the central formulas. Our development of the perturbation theory differs somewhat from Peliti's. Most of the applications discussed are also new.
The balance of this article is structured as follows. In Sec. II we derive the path-integral representation, starting from the master equation. Sec. III presents an application to the simple decay process, and expressions for two-time joint probabilities. In Sec. IV we begin our discussion of diagrammatic perturbation theory for the probability generating function, which is illustrated with a pedagogical example. This is extended in Sec. V where we analyze the birth-and-death process using perturbation theory. In Sec. VI a perturbation expansion for moments of the distribution is developed, which turns out to be much simpler than that for the full generating function. This method is applied to the Malthus-Verhulst process in Sec. VII. In Sec. VIII we illustrate another application of the formalism, showing how the path-integral description for a lattice of coupled Malthus-Verhulst processes leads, in the continuum limit, to a field theory for directed percolation. Sec. IX presents a brief summary.
2 From the master equation to a path integral
In this section we recapitulate Peliti's derivation of the path integral mapping. We consider Markov processes in continuous time, and with a discrete state space n = 0,1,2,.... (We may think of n as the size of a certain population.) The probability pn(t) of state n at time t is governed by the master equation [10,11,12]:
where wmn is the rate for transitions from n to m. (We study stationary stationary Markov processes, i.e., time-independent transition rates.)
We now associate a vector |nñ in a Hilbert space with each state n, and for convenience define the inner product so:
The identity may then be written as
(All sums run from zero to infinity unless otherwise specified.) A probability distribution fn (fn> 0, ånfn = 1), may be represented as a linear combination of basis states:
The Hilbert space formalism is useful because it provides a simple way to express the evolution in terms of creation (p) and annihilation (a) operators, which we define via:
These relations imply that [a,p] = 1. (Note that while we make use of many pieces of notation familiar from quantum mechanics, there are fundamental differences. Expected values, for example, are linear, not bilinear, in |fñ.) The mean population size is given by
where
is the the projection onto all possible states.
Central to our analysis will be the probability generating function (PGF),
We denote the PGF corresponding to state |fñ as f(z) = ånfnzn. (Note that f(1) = 1 by normalization.)
Next consider the inner product between states |fñ and |yñ:
We can write this in terms of the corresponding PGFs if we note the identity
which is readily proved, integrating by parts. (Unless otherwise specified, all integrals are over the real axis.) Then we have
where we used the integral representation of the d function.
For birth-and-death processes, it is always possible to write the master equation in terms of an evolution operator L composed of creation and annihilation operators (specific examples are considered below). The master equation then takes the form
and has the formal solution
This evolution has its analog in the space of probability generating functions; it is for the analog of the operator Ut in the PGF representation that we shall develop a path-integral expression. To do this we define, for any operator A in the Hilbert space, a function called its kernel:
where Am,n = ám|A|nñ are the matrix elements of A. Suppose |yñ = A |fñ. The PGF corresponding to |yñ is
We shall also require an expression for the kernel of a product of a pair of operators, A and B:
Given an operator A, we may put it in normal order by commuting all creation operators to the left of all annihilation operators; it will then have the form:
With this we associate the normal kernel:
The ordinary and normal kernels are related via:
[To prove Eq.(19) we show that the coefficients of zmzn on the right- and left-hand sides are equal. The coefficient on the r.h.s. is:
while on the l.h.s. it is simply Am,n/(m!n!). The matrix element, however, may be written so:
which establishes the identity.]
We may now develop a path-integral representation for Ut(z,z). To begin, recall Trotter's formula, which allows us to write:
Each factor in the product has a corresponding kernel, which, using Eq. (19), can be written
with (z,z) the normal kernel of the evolution operator. Now using Eq. (16), we have
or, rearranging,
Finally, we let N ® ¥ with t¢ = (k/N)t, hk® y(t¢) and hk¢® y¢(t¢) , to obtain a path-integral expression for Ut (z,z):
where the dot denotes a time derivative. The functional integrals over y(s) and y¢(s) are for 0 < s < t, with boundary conditions y(0) = z and iy¢(t) = z; y and y¢ are real. The symbol òDy is defined by the limiting process in Eq. (23). [The reader may wonder at this point what has become of the factors of 2 p in the denominator of Eq. (23). The answer is that the prefactor, which for the moment is undefined, will be fixed via normalization.] Note also that the first term in the argument of the exponential could be written more precisely as iy¢(t¢) (t¢), i.e., the time derivative is evaluated at t¢- with ® 0 from above. While the function y¢(t¢) has no obvious physical significance, we shall see that y is closely related to the random variable n(t) in the birth-and-death process.
The kernel Ut (z,z) has two principal uses. First, from Eqs. (15) and (16), we see that Ut provides a mapping between PGFs at different times:
We make considerable use of this relation in the examples that follow. Evaluating the integral is particularly simple if the initial distribution is Poisson. Then F0(z) = ep(z1), and
Another simple case is pn(0) = dn,n0 (exactly n0 individuals initially), corresponding to F0(z) = zn0, which yields
Setting z = 0 in the above expressions, we obtain the probability of the state n = 0, which in many cases is absorbing, so that the survival probability is given by Ps (t) = 1 p0 (t).
Independent of its probability interpretation, the kernel has a second utility: we may treat the argument of the exponential as an effective action. In the continuum limit of a system with many degrees of freedom, the function y(t) becomes a classical field, y(x,t) (similarly for y¢), leading to a field theory for a Markov process originally described by transition rates for particles on a lattice. (An example is discussed in Sec. VIII.) With such a field theory in hand, we can apply methods such as the renormalization group to study critical behavior. The effective action is known once we construct the evolution operator L; at no point do we need to write a Langevin equation or stipulate noise properties.
3 Decay Process
As a simple example we consider exponential decay, i.e., the Markov process with transition rates wm,n = wn dm,n1. The evolution operator is
Since this is in normal order, we have
and thus,
It is convenient to transform away the linear term in the action via the change of variable i = iy¢1. The first term in the integral then contributes the additional (boundary) terms y(t) + y(0), and using y(0) = z we have
This can be evaluated exactly. Recalling that ò dweiwz = 2pd(z) , we see that the functional integral
imposes the condition
for 0 < s < t, and so y(t) = y(0) ewt = zewt. In other words the functional integral over yields a product of d-functions:
which when inserted in Eq. (31) gives
where C represents the as yet undetermined normalization factor. Using this in Eq. (25) results in
where in the final line we set C = 1 to satisfy the normalization condition, Ft(1) = 1. If there are exactly n particles at time zero, F0(z) = zn, and
which on expanding yields
as expected. For a Poisson initial distribution we find
corresponding to a Poisson distribution whose mean decays exponentially: p(t) = pewt.
3.1 Joint probabilities
It is useful to extend the formalism to joint probabilities, i.e., for the values of the process at different times. For t1> t2, let P(n1,t1;n2,t2|n0,0) be the probability of state n1 at time t1and n2 at t2, given n0 at time 0. The generating function for the joint probability is
where in the second line we used the Markov property and in the third, stationarity.
Let us modify slightly our notation for the one-time generating function, to include the initial condition, P(n,0) = fn (0) = dn,n0:
Then we have
where, in the second line, we used Eq. (25). Now, noting that F(iz¢,0|n2) = (iz¢)n2, we have
Thus we have a formula analogous to Eq. (25), for the two-time generating function. The generalization to n times is straightforward.
For the decay process, Eq. (43) reads:
which, after integrations over d-functions, yields:
From this we readily obtain the expectation:
and
The covariance is then
For
t
1 =
t
2 we find áá
n(
t) ññ =
n
0
e
wt [1
e
wt] , which is the variance of a binomial random variable [see Eq. (38)].
4 Perturbation Theory
The preceding example illustrates the general procedure for calculating probabilities, but no perturbative treatment was needed as the action was purely bilinear in the fields. An example of perturbation theory is afforded by analyzing the Malthus-Verhulst process, in which each individual has a rate to reproduce, and a rate of m+ n(n1) to die, if the total population is n. (The term µ n represents saturation, so that the population does not grow without limit even if > m.) For this process the kernel is
since the corresponding evolution operator has the property
which evidently reproduces the rates defining the process.
Using Eq. (25) we have (with º iy¢),
As before, we eliminate the linear term by changing variables, = iy¢1, yielding
where w = m . Separating the bilinear part of the action from the cubic and quartic terms, we have
where
Evidently SI must be treated perturbatively. It is interesting to note that while the process with n = 0 admits an exact solution, the cubic term generates a perturbation series in the present formalism. When we expand eSI we obtain a series of functional integrals over y and of various products of these fields times the "Gaussian" factor eS0. Consider the basic contraction:
Let
where in the second line we integrated by parts and used the boundary conditions y(0) = z and (t) = z1. Since the contraction [y(t1) (t2)] is the functional integral of y(t1) (t2) , it is convenient to introduce the operator via the relation
The fact that
suggests that take the form
where is a boundary term. Using Eq. (58) and integrating by parts yields
Consistency with Eq. (57) then requires that
which has the 'causal' solution
leading to
as may be verified directly.
To evaluate
we make use of our earlier result, Eq. (33), to evaluate the integral over :
Now perform a functional integration by parts so that d/dy(t¢¢) operates on y(t1) . This yields
We may now integrate over y to obtain
where
In light of the discussion following Eq. (24), we see that in Eq. (57) and the subsequent expressions, t2 should be interpreted as t, which means that in the expectation, Eq. (67), we should take Q(0) º 0.
To see what our results mean, consider the very simple situation of
I = uy, which is equivalent to exponential decay with a rate w¢ = w + u. Then,
and
consistent with the exact result, Ut(z,z) = ez[1 + (z1)ewteut]. (Note that obtaining the correct result depends on using Q(0) = 0.)
The preceding analysis exposes a curious feature of this formalism. Normally when we perturb about a Gaussian model, the expectation of the field áyñ = 0. Here, by contrast, we have
while from Eq. (65) we have
This suggests that it would simplify matters if we were to introduce new fields:
and
which by construction have expectation zero. The boundary conditions on y and imply that j(0) = (t) = 0.
Noting that ¶t¢ + w annihilates ewt¢, we have
where in the second line we integrated by parts and used the fact that ¶t¢ + w annihilates ewt¢. Integrating by parts once again gives
which, when inserted in Eq. (53), yields
The change of variable, then, yields the unperturbed solution, U, automatically, and eliminates the boundary term from the exponential. The normalization is such that the functional integral is unity if SI = 0.
We may now repeat the analysis of the basic contraction by letting
and defining such that
One readily verifies that
which leads directly to
Then an average of n fields (with equal numbers j's and 's) may be written
with áj(t1) (t2) ¼j(tn1) (tn)ñ given by the sum of all distinct products of pairwise contractions of the form of Eq. (81).
Let us analyze the simple quadratic "perturbation" using the "j" representation. Introducing the new fields we have
which gives us
where
and
Let Ut(z,z) = U¢ (z,z) (z,z). The latter factor has the expansion:
where jjº j(tj), etc. The n = 0 term is unity (by normalization), while the n = 1 term vanishes, since áj(t¢)ñ = á(t¢)ñ = áj(t¢) (t¢)ñ = 0.
For n = 2, of the nine possible terms in , only one survives the functional integral, giving
Next consider the n = 3 contribution. Once again there is but a single nonzero term, since we must choose b1j1 for the first factor (there is no j with time > t1, with which to contract 1), c3
3 for the final factor (there is no with time < t3, i.e., no way to contract j3), while the middle factor must be 2j2. Writing this contribution in terms of functional derivatives we have
Upon performing the two functional integrations by parts, d/dj(t¢¢) may act on j1 and d/dj(t¢¢¢) on j2, or vice-versa. In the former case we obtain the factor
The d-functions can both be satisfied since t1Î [t2,t] and t2Î [t3,t]. On the other hand, the second alternative gives the factor
since t¢¢ = t2 does not fall within the range of integration. Combining this result with the other factors the third-order term becomes z(z1) ewt (ut)3/3!.
We now introduce a diagrammatic notation which will prove essential in this and subsequent analyses. The three terms in
I are represented by vertices, so:
A contraction between ji and j (with ti > tj due to the factor Q(ti tj)) is represented by joining the line exiting vertex j with that entering vertex i. Thus the n = 2 and n = 3 terms correspond to the diagrams:
respectively. Associated with the line connecting vertices i and j (with ti > tj), is the factor ew(titj).
Each diagram constructed according to the following rules corresponds to a term in the expansion of .
1) Draw m > 1 vertices of type ·-¬ and an equal number of type -¬·. Add p > 0 vertices of type -¬·-¬, for a total of 2m+p vertices.
2) Form all possible unlabelled diagrams, i.e., distinct connections in which each outgoing line is joined with an ingoing line.
3) Generate all distinct labellings of each unlabelled diagram by assigning an index (time) to each vertex, such that the arrows always point from the larger to the smaller index.
4) A diagram consists of one or more connected parts; its contribution to is the product of factors associated with these connected parts.
In the example under consideration, the factor associated with the jth connected part (having vj vertices) is
We now make use of a well known theorem in diagrammatic analysis (see, e.g., £ 8.3 of Ref. [14]), which in the present context states that ln
t is given by the sum of all contributions associated with connected diagrams only. In the present case there is exactly one connected diagram at each order n > 2, yielding
which, when inserted in Eq. (84), gives the exact result,
While the expansion of ln
t in terms of connected graphs represents a considerable simplification, it is possible, in principle, to evaluate t directly, without transforming to the fields j and . Consider, for example, the second order term in the expansion of t , Eq. (70):
Noting that the functional derivative w.r.t. y() gives zero, the above expression is seen to be
which is precisely the (u2) contribution in the expansion of t, Eq. (94). Clearly, this analysis is considerably more cumbersome than the connected diagram expansion.
5 Birth-and-death process
Next we consider the more complex, but still exactly soluble problem of a birth-and-death process without saturation, i.e., Eq. (54) with n = 0. Our first step, as before, is to rewrite the perturbation in terms of the variables j and :
In the second line, the final term is independent of j and . Integrating it, and combining it with the usual prefactor U, we have for this problem
The remaining terms will again be analyzed using diagrams. There are five vertices, as shown in Fig. 1 (note that the three we have encountered before carry different factors here).
The first vertex in Fig. 1 ("terminal") must bear the lowest index of a given branch, since no lines exit from it. The second and fifth can only appear in the interior of a diagram, while the remaining two can only appear as the nth of an n-vertex diagram, since no lines enter. Thus the expansion is not as complicated as it might at first appear. The lowest-order diagram is again which now corresponds to
At third and higher orders there are graphs with bifurcations, the simplest being the one shown in Fig. 2, corresponding to
The factor of 2 arises because there are two ways to contract the lines exiting vertex 3. Thus all the vertices in this problem, except for the terminal b2j, carry a factor of 2, either explicitly, or due to the combinatorial factor. The other third-order connected diagram is
The two diagrams differ only by their numerical prefactors, the sum of which just cancels the 1/3! that comes from the integrations. We show in the AppendixAppendix), and so ), and a factor We aim to show that the sum of all numerical weights associated with n-vertex diagrams in the expansion of lnt for the birth-and-death process is n!. To see this, consider an arbitrary n,b-diagram, i.e., one having n vertices, b of them bifurcations. There is a single factor of z, associated with vertex n. Such a diagram will have b+1 terminal vertices (carrying factors b), and so n2b1 nonterminal vertices that are not bifurcations, which carry factors bj. Each bk carries a factor of ewt (z1) for a total of n+1 such. Next consider the factors ewti for each vertex i. If i is terminal, there are two such factors (from b), and a factor ewti due to the line entering the vertex. If i < n is neither terminal nor a bifurcation, there is a factor ewti from bi, while the exponential factors associated with the lines entering and leaving the vertex cancel. If i < n is a bifurcation, there is no factor bi, but there is again a net factor of ewti as two lines exit, while only one enters vertex i. Consider vertex n. If it is a bifurcation, then the net factor is ewtn (two lines exiting, factor ewtn in cn). The same holds if vertex n is not a bifurcation (single line exiting, factors in bi and ci cancel). In summary, there is a factor ewti associated with each vertex. Finally, all vertices carry a factor of 2, except for the terminal ones, leading to an overall factor of 2nb1. Combining all of these observations, the contribution due to a given labelled n,bdiagram is that this pattern continues at higher orders, i.e., that the sum of all numerical factors due to n-vertex diagrams is n!, so that
Combining this with Eq. (98), we obtain the exact result [13],
6 Expansion of moments
While the transformation to variables j and , with expectation zero, yields a certain simplifcation, it also causes a proliferation in the number of vertices, so that in many cases it may be advantageous to work with the original variables y and . We have seen that in this case the expansion of Ut becomes quite complicated. Fortunately, it is still possible to derive relatively simple expressions for the moments ánr (t)ñ of the distribution.
Note that from the definition of Ft(z) we have
and that, in general, the r-th factorial moment, ánrñfº án(n1)(n2) ¼(nr+1) ñ, is given by:
Using Eq. (25) we then have
If Ut is of the form of Eq. (53), then
(Note that the prefactor ez is just U0 for z = 1.) Thus the expectation of the variable y is closely related to that of n(t) itself.
In case n(0) = n0 (so that F0(z) = zn0), we have, using Eq. (27),
Evaluation of the above expression is facilitated by noting that
If the initial distribution is Poisson with parameter p, then Eq. (26) implies that
It is instructive to evaluate the above expressions first for the simple decay process, SI = 0. In this case
and a short calculation shows that for the initial condition n(0) = n0,
and
as expected. For the Poisson initial distribution, one finds that
consistent with the probability distribution pn(t) being Poisson with parameter pewt.
Next we consider the moments for the birth-and-death process. In this case
corresponding to a vertex that branches to the left. To evaluate U (z), we expand eSI, generating diagrams consisting of vertices and a single "sink", with r lines entering, to the left of all vertices. Each line entering a node (i.e., a vertex or the sink), corresponds to a variable y(t¢); if not contracted, it contributes a factor of zewt¢. (We call such uncontracted lines "external". Recall that for the birth-and-death process, w º m .) On the other hand, lines exiting a vertex correspond to variables, and must be contracted with a y line, because [] = 0 when z = 1. With z = 1, the basic contraction is:
Thus the perturbation series for a given moment is much simpler than the full series for Ut (z,z).
Consider the case r = 1. The sink has only a single line entering, so diagrams with n > 1 vertices all give zero, and the only nonzero contribution to (¶Ut /¶z)z = 1 is for n = 0 (the sink itself), which is zez ewt.
For r = 2 the sink has two lines. There is then the n = 0 contribution, z2 ez e2wt, and the 1vertex diagram shown in Fig. 3. Including the combinatorial factor of 2, its contribution is
Using these results, we readily find that
For n(0) = n0, the variance is
for w ¹0, while for w = 0 we have
For a Poisson initial distribution,
which becomes p(1 + 2 t) in case w = 0. We see that at long times, the variance decays (grows) exponentially, for w > 0 (w < 0), and that for w = 0 the variance grows ~ t, reflecting the diffusive character of the process in this case. Higher-order moments may be evaluated similarly; only diagrams with n < r 1 vertices contribute to ánrñf.
7 Malthus-Verhulst Process
We now return to the Malthus-Verhulst process; the three vertices associated with SI [see Eq. (54)] are shown in Fig. 4. The series for U (z) begins with the "sink" term, zr ez erwt, and then includes diagrams with n > 1 vertices. (Note that the n-th order contributions carry factors of (n)s
ns, with s variable.) Based on our diagrammatic analysis of the birth-and-death process, we can formulate the following simple rules for evaluating the nth order contribution to U:1) Draw all diagrams consisting of n vertices and a single r-line sink to the left of all vertices. Each line exiting a vertex must be contracted with a line entering another vertex to the left.
2) The sink has time variable t; assign time variables t1,...,tn to the vertices from left (nearest the sink) to right.
3) For each external line include a factor of zew, where is the time variable of the associated vertex (if the external line is attached to the sink, = t). For each internal line from vertex j to vertex i, (i < j) include a factor of ew(ti tj).
4) Include the factors ( or n) associated with each vertex, and the combinatorial factor associated with the number of ways of realizing the contractions. Finally, integrate over the time variables t1,...,tn, with
(Note that by fixing the time ordering we effectively include the factor 1/n! that comes from expanding eSI.)
Let us consider some examples in the evaluation of ánñ; the low-order diagrams are shown in Fig. 5. For n = 1 we have only diagram (a), whose contribution is:
At second order we have the following contributions. Diagram (b):
Diagram (c):
Diagram (d):
For fixed n(0) this yields the expansion:
For r = 2 a similar calculation results in
The evaluation of diagrams can be simplified somewhat by considering the Laplace transform of ánr ñf:
We have seen that each line entering a node with time variable carries a factor of ew, and each line exiting the node carries a factor of ewt. Letting i be the number of lines entering node i less the number exiting, we can write the time integration factor for the general diagram so:
Inverting the order of the integrations, the Laplace transform of the above expression becomes:
7.1 A numerical example
The diagramatic expansion for the Malthus-Verhulst process can be extended, and simplified (by identifying the so-called reducible diagrams [14]), but such analysis is beyond the scope of this article and will be defered to a future work. We close this section with an application of the second-order expansion for the mean population, Eq. (117). Note that for short times, the factor multiplying n0ewt may be written as
with A = n(n0 1) and B = n[(n0 1)(n0 3)n ] . Evidently, our analysis generates an expansion in powers of t, whose convergence would seem to require that nn0 t << 1.
If we wish to extend the range of validity of the expansion, we must transform to a new variable that remains finite as t ® ¥; a glance at Eq. (117) suggests that we use y = 1 ewt as the new variable. (This transformation is very useful in analysis of series, for example in the study of random sequential adsorption [15].) In fact, f(t) is readily expressed in terms of y:
with = A/w and = B/w2 To proceed, we form a Padé approximant [16,17] to the power series in y:
Equating coefficients of y and of y2, we find
and
(Here we must note that a Padé approximant to a series of three terms can only serve as a very rough approximation!)
We apply this expression to the Malthus-Verhulst process with parameters m = 1, = 0.5, and n = 0.1. In Fig. 6 we compare the numerical solution of the master equation with the perturbation theory expression, án(t) ñ = n0 ewtf(t) (with f represented by the Padé approximant), and with the simple exponential decay, n0 ewt, for n0 = 3 and n0 = 10. The perturbation expansion, using the Padé expression, yields a good approximation to the correct value, despite the very short series used, and the fact that for n0 = 10 we have nn0 = 1, so that the time series has only a small radius of convergence.
8 Coupled Malthus-Verhulst processes
In this section we show how the path integral formalism can be applied to processes on a lattice, leading to a field theory in the continuum limit. The process of interest is a lattice of coupled Malthus-Verhulst processes, with diffusive exchange of particles between neighboring cells. This process describes the dynamics of a population distributed in space, with the individuals performing random walks, allowing the population to spread. In the continuum limit, the process is of great interest as a field theory for directed percolation [6].
At each site r of the lattice, there is a process nr(t) whose evolution is given by Eq. (50). The new feature in this model is diffusion, represented by the operator:
where D is the diffusion rate and åe is over the vectors from a given site to its nearest neighbors.
The evolution kernel Ut is a function of the variables zr and zr defined at each site; we use {z} and {z} to denote these sets of variables. To write Ut for this process we add the diffusive contribution to the Malthus-Verhulst part, Eq. (51), to find
where Dyr = åe (yr+e yr) is the discrete Laplacian, and the functional integrals are now understood to include the funcions yr and r at each site. Eliminating the linear term as usual, by letting = iy¢1, and introducing w = m , we obtain
Up to this point, our expression for Ut is exact. We now make a number of simplifications, leading to an effective action for the lattice of coupled Malthus-Verhulst processes. (Each deserves a careful justification, but we shall not enter into such questions here.)
i) We drop the boundary term, which should not influence stationary properties.
ii) We take the continuum limit, so that yr(t)® y(x, t), and Dyr ® Ñ2 y.
iii) We discard the term µ
2 y2, which turns out to be irrelevant to the scaling behavior near the critical point [5].Under these approximations the argument of the exponential in Eq. (123) becomes the effective action
This is the action corresponding to directed percolation [6] or Reggeon field theory (a particle physics model with the same formal structure as the continuum theory of a spatially extended population) [18,19]. This theory has been analyzed using renormalization group techniques, to show that the upper critical dimension dc = 4, and to derive expressions for critical exponents in an expansion in = 4 d [20,21]. While such developments lie beyond the scope of this article, we can get a qualitative understanding of the physics represented by S by ignoring, for the moment, the term µ y2 . Functional integration over then imposes the following partial differential equation as a constraint on y:
This is the mean-field theory of directed percolation and allied models [22]. We see that this equation correctly predicts an absorbing state, y = 0. For w > 0 the solution flows to this state, regardless of the initial condition. For w < 0, however, another stationary state appears: the uniform solution y = |w|/n. Thus w = 0 marks a continuous phase transition. The field y(x,t) may be interpreted as a local population density. This mean-field description is readily shown to yield the critical exponents b = 1 (for the order parameter áyñ), n|| = 1 (for the relaxation time), and n^ = 1/2 (for the correlation length) [22].
In this case, the effect of the neglected term in the action, y2 , can be represented by a noise term h(x,t) in Eq. (125), which now becomes a Langevin equation or stochastic partial differential equation:
where h is a Gaussian noise with áh(x,t) ñ = 0 and autocorrelation proportional to the local density [6]:
(note that in this way the noise respects the absorbing state). In the presence of noise, the critical point is renormalized from its mean-field value of wc = 0, and, more significantly, the critical exponents take non-mean-field values for d < 4. As we have noted, the transcription of a stochastic model to a Langevin equation is not always straightforward. For this reason, the path-integral formalism discussed in this article is especially valuable: it allows one to construct an action (that is, the starting point for a renormalization group analysis), without having to postulate noise properties.
9 Summary
We have reviewed the formalism, based on the work of Peliti, that maps a Markov process to a path integral representation, and have presented detailed examples of its application to birth-and-death processes. We show, in particular, how the exact solutions for the decay and simple birth-and-death process can be recovered, and derive a perturbation expansion for moments in the Malthus-Verhulst process. Finally, we show how the evolution kernel for a lattice of coupled Malthus-Verhulst processes leads, in the continuum limit, to a field theory for directed percolation. As the study of nonequilibrium processes grows, we expect the methods discussed here to attract ever-greater interest.
Acknowledgments
We are grateful to Migual A. Muñoz and Adélcio Carlos de Oliveira for helpful comments. This work was supported by CNPq.
To find the contribution µ
n, it remains to evaluate the sum of all such terms, that is, to determine
where the sum is over all distinct labelled diagrams of n vertices. A labelled n,bdiagram has a weight of 2nb1. Let W(n,b) be the sum the weights of all n,bdiagrams. Define the degree of a vertex as the number of lines that exit from it, so that a terminal vertex has degree 0, bifurcations have degree 2, and all other vertices have degree 1. (In this model there are no vertices with degree > 2.)
Now, given a labelled n,bdiagram, we can generate a set of distinct n+1vertex labelled diagrams by the following recipe:
i) Relabel the vertices 1,...,n as 2,...,n+1.
ii) Attach a new vertex ('1') to any vertex of degree less than 2.
It is easy to see that (1) each choice for attaching the new vertex generates a different diagram; (2) the sets generated by different nvertex labelled diagrams are mutually disjoint; (3) applied to the complete set of nvertex labelled diagrams, the procedure generates the complete set of n+1vertex diagrams. When we attach the new vertex to one of degree 0, the number of bifurcations does not change, so the new diagram has an additional factor of 2 in its weight. If on the other hand we attatch the new vertex to one of degree 1, we generate a new bifurcation, and the weight remains unchanged. Recalling that an n,bdiagram has b+1 vertices of degree zero, and n2b1 of degree 1, we have for n > 2 the following recurrence relation:
Starting with W(2,0) = 2 (and W(2,j) = 0 for j > 0), we readily find W(3,0) = 4, W(3,1) = 2, and then W(4,0) = 8, W(4,1) = 16, and so on.
To solve for W(n) in general, we introduce a generating function
(Since W(n) appears to grow like n! we need the factorial for g to have a nonzero radius of convergence. Without it, we get a function with an essential singularity at x = 0, as the reader may verify!) For purposes of analysis, it is convenient to define W(1,0) º 1 and W(1,j) º 0 for j > 0, which is completely consistent with the recurrence relation for n = 2. (There is, of course, no diagram with n = 1.) For n = 0, W vanishes identically, so we must add a source term in Eq. (130). The modified recurrence relation, valid for n = 1, 2, 3,... and b > 0, is
Multiply this relation by xn yb/n! and sum over n and b. Letting m = n1 and rearranging, one finds
Now let b¢ = b1 in the terms multiplying W(m,b1). After a simple rearrangement (and dropping the primes) we have
We require
For y = 1, Eq. (134) immediately yields G(x) = x/(1x), implying W(n) = n!.
References
[1] P. C. Martin, E. D. Siggia, and H. A. Rose, Phys. Rev. A8, 423 (1978).
[2] C. De Dominicis, J. Physique 37, 247 (1976).
[3] M. Doi, J. Phys. A9, 1465, 1479 (1976).
[4] P. Grassberger and M Scheunert, Fortschr. Phys. 28, 547 (1980).
[5] H. K. Janssen, Z. Phys. B23, 377 (1976).
[6] J. L. Cardy, Scaling and Renormalization in Statistical Physics, (Cambridge University Press, Cambridge 1996), Ch. 10.
[7] C. W. Gardiner and S Chaturvedi, J. Stat. Phys. 17, 429 (1977); 18, 501 (1978); D. Elderfield and D. D. Vvedensky, J. Phys. A18, 2591 (1985); M. Droz and A. McKane, J. Phys. A27, L467 (1994).
[8] L. Peliti, J. Physique 46, 1469 (1985).
[9] M. J. Howard and U. C. Täuber, J. Phys. A30, 7721 (1997).
[10] N. G. van Kampen, Stochastic Processes in Physics and Chemistry, (North-Holland, Amsterdam, 1992).
[11] C. W. Gardiner, Handbook of Stochastic Methods, (Springer-Verlag, Berlin, 1990).
[12] T. Tomé and M. J. de Oliveira, Dinâmica Estocástica e Irreversibilidade, (EdUSP, São Paulo, 2001).
[13] J. Marro and R. Dickman Nonequilibrium Phase Transitions in Lattice Models (Cambridge University Press, Cambridge, 1999), Ch. 1.
[14] J. J. Binney, N. J. Dowrick, A. J. Fisher, and M. E. J. Newman, The Theory of Critical Phenomena, (Oxford University Press, Oxford, 1992).
[15] R. Dickman, I. Jensen, and J.-S. Wang, J. Chem. Phys. 94, 8252 (1991).
[16] G. A. Baker, Jr. Essentials of Padé Approximants, (Academic PRess, New York, 1975).
[17] A. J. Guttmann, in Phase Transitions and Critical Phenomena, Vol. 13, C. Domb and J. Lebowitz, eds. (Academic Press, New York, 1989).
[18] J. L. Cardy and R. L. Sugar, J. Phys. A13, L423 (1980).
[19] P. Grassberger and K. Sundremeyer, Phys. Lett. B77 220, (1978); P. Grassberger and A. de la Torre, Ann. Phys. (NY) 122, 373 (1979).
[20] H. K. Janssen, Z. Phys. B42, 151 (1981).
[21] F. van Wijland, K. Oerding, and H. J. Hilhorst, Physica A251, 179 (1998).
[22] See Chapter 6 of Ref. [13].
Received on 27 August, 2002
- [1] P. C. Martin, E. D. Siggia, and H. A. Rose, Phys. Rev. A8, 423 (1978).
- [2] C. De Dominicis, J. Physique 37, 247 (1976).
- [3] M. Doi, J. Phys. A9, 1465, 1479 (1976).
- [4] P. Grassberger and M Scheunert, Fortschr. Phys. 28, 547 (1980).
- [5] H. K. Janssen, Z. Phys. B23, 377 (1976).
- [6] J. L. Cardy, Scaling and Renormalization in Statistical Physics, (Cambridge University Press, Cambridge 1996), Ch. 10.
- [7] C. W. Gardiner and S Chaturvedi, J. Stat. Phys. 17, 429 (1977); 18, 501 (1978);
- D. Elderfield and D. D. Vvedensky, J. Phys. A18, 2591 (1985);
- M. Droz and A. McKane, J. Phys. A27, L467 (1994).
- [8] L. Peliti, J. Physique 46, 1469 (1985).
- [9] M. J. Howard and U. C. Täuber, J. Phys. A30, 7721 (1997).
- [10] N. G. van Kampen, Stochastic Processes in Physics and Chemistry, (North-Holland, Amsterdam, 1992).
- [11] C. W. Gardiner, Handbook of Stochastic Methods, (Springer-Verlag, Berlin, 1990).
- [12] T. Tomé and M. J. de Oliveira, Dinâmica Estocástica e Irreversibilidade, (EdUSP, Săo Paulo, 2001).
- [13] J. Marro and R. Dickman Nonequilibrium Phase Transitions in Lattice Models (Cambridge University Press, Cambridge, 1999), Ch. 1.
- [14] J. J. Binney, N. J. Dowrick, A. J. Fisher, and M. E. J. Newman, The Theory of Critical Phenomena, (Oxford University Press, Oxford, 1992).
- [15] R. Dickman, I. Jensen, and J.-S. Wang, J. Chem. Phys. 94, 8252 (1991).
- [16] G. A. Baker, Jr. Essentials of Padé Approximants, (Academic PRess, New York, 1975).
- [17] A. J. Guttmann, in Phase Transitions and Critical Phenomena, Vol. 13, C. Domb and J. Lebowitz, eds. (Academic Press, New York, 1989).
- [18] J. L. Cardy and R. L. Sugar, J. Phys. A13, L423 (1980).
- [19] P. Grassberger and K. Sundremeyer, Phys. Lett. B77 220, (1978);
- P. Grassberger and A. de la Torre, Ann. Phys. (NY) 122, 373 (1979).
- [20] H. K. Janssen, Z. Phys. B42, 151 (1981).
- [21] F. van Wijland, K. Oerding, and H. J. Hilhorst, Physica A251, 179 (1998).
Appendix
), and so ), and a factorWe aim to show that the sum of all numerical weights associated with n-vertex diagrams in the expansion of lnt for the birth-and-death process is n!. To see this, consider an arbitrary n,b-diagram, i.e., one having n vertices, b of them bifurcations. There is a single factor of z, associated with vertex n. Such a diagram will have b+1 terminal vertices (carrying factors b), and so n2b1 nonterminal vertices that are not bifurcations, which carry factors bj. Each bk carries a factor of ewt (z1) for a total of n+1 such. Next consider the factors ewti for each vertex i. If i is terminal, there are two such factors (from b), and a factor ewti due to the line entering the vertex. If i < n is neither terminal nor a bifurcation, there is a factor ewti from bi, while the exponential factors associated with the lines entering and leaving the vertex cancel. If i < n is a bifurcation, there is no factor bi, but there is again a net factor of ewti as two lines exit, while only one enters vertex i. Consider vertex n. If it is a bifurcation, then the net factor is ewtn (two lines exiting, factor ewtn in cn). The same holds if vertex n is not a bifurcation (single line exiting, factors in bi and ci cancel). In summary, there is a factor ewti associated with each vertex. Finally, all vertices carry a factor of 2, except for the terminal ones, leading to an overall factor of 2nb1. Combining all of these observations, the contribution due to a given labelled n,bdiagram is
Publication Dates
-
Publication in this collection
23 Apr 2003 -
Date of issue
Mar 2003
History
-
Received
17 Aug 2002