Path integrals and perturbation theory for stochastic processes

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}

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 p_n(t) of state n at time t is governed by the master equation [10,11,12]:

where w_mn 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 f_n (f_n> 0, å_nf_n = 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) = å_nf_nzⁿ. (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 U_t 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 A_m,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 z^mzⁿ 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 A_m,n/(m!n!). The matrix element, however, may be written so:

which establishes the identity.]

We may now develop a path-integral representation for U_t(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, h_k® y(t¢) and h_k¢® y¢(t¢) , to obtain a path-integral expression for U_t (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 U_t (z,z) has two principal uses. First, from Eqs. (15) and (16), we see that U_t 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 F₀(z) = e^p(z–1), and

Another simple case is p_n(0) = d_n,n0 (exactly n₀ individuals initially), corresponding to F₀(z) = z^n₀, 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 P_s (t) = 1 – p₀ (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 w_m,n = wn d_m,n–1. 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 ò dwe^iwz = 2pd(z) , we see that the functional integral

imposes the condition

for 0 < s < t, and so y(t) = y(0) e^–wt = ze^–wt. 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, F_t(1) = 1. If there are exactly n particles at time zero, F₀(z) = zⁿ, 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) = pe^–wt.

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 t₁> t₂, let P(n₁,t₁;n₂,t₂|n₀,0) be the probability of state n₁ at time t₁and n₂ at t₂, given n₀ 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) = f_n (0) = d_n,n0:

Then we have

where, in the second line, we used Eq. (25). Now, noting that F(iz¢,0|n₂) = (iz¢)^n₂, 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

₁ =

t

₂ we find áá

n(

t) ññ =

n

₀

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(n–1) 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 S_I 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 e^–S_I we obtain a series of functional integrals over y and of various products of these fields times the "Gaussian" factor e^–S₀. Consider the basic contraction:

Let

where in the second line we integrated by parts and used the boundary conditions y(0) = z and (t) = z–1. Since the contraction [y(t₁) (t₂)] is the functional integral of y(t₁) (t₂) , 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(t₁) . 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, t₂ 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, U_t(z,z) = e^{z[1 + (z–1)e–wte–ut]}. (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 e^–wt¢, we have

where in the second line we integrated by parts and used the fact that –¶_t¢ + w annihilates e^wt¢. 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 S_I = 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(t₁) (t₂) ¼j(t_n–1) (t_n)ñ 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 U_t(z,z) = U¢ (z,z) (z,z). The latter factor has the expansion:

where j_jº j(t_j), 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 b₁j₁ for the first factor (there is no j with time > t₁, with which to contract ₁), c₃

₃ for the final factor (there is no with time < t₃, i.e., no way to contract j₃), while the middle factor must be ₂j₂. Writing this contribution in terms of functional derivatives we have

Upon performing the two functional integrations by parts, d/dj(t¢¢) may act on j₁ and d/dj(t¢¢¢) on j₂, or vice-versa. In the former case we obtain the factor

The d-functions can both be satisfied since t₁Î [t₂,t] and t₂Î [t₃,t]. On the other hand, the second alternative gives the factor

since t¢¢ = t₂ does not fall within the range of integration. Combining this result with the other factors the third-order term becomes z(z–1) e^–wt (–ut)³/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 j_i and _j (with t_i > t_j due to the factor Q(t_i –t_j)) 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 t_i > t_j), is the factor e^{–w(t_i–t_j)}.

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 j–th connected part (having v_j 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 (u²) 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 b²j, 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 ^Appendix 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 U_t becomes quite complicated. Fortunately, it is still possible to derive relatively simple expressions for the moments án^r (t)ñ of the distribution.

Note that from the definition of F_t(z) we have

and that, in general, the r-th factorial moment, án^rñ_fº án(n–1)(n–2) ¼(n–r+1) ñ, is given by:

Using Eq. (25) we then have

If U_t is of the form of Eq. (53), then

(Note that the prefactor e^z is just U⁰ for z = 1.) Thus the expectation of the variable ^y is closely related to that of n(t) itself.

In case n(0) = n₀ (so that ^F₀(z) = z^n₀), 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, S_I = 0. In this case

and a short calculation shows that for the initial condition n(0) = n₀,

and

as expected. For the Poisson initial distribution, one finds that

consistent with the probability distribution p_n(t) being Poisson with parameter pe^–wt.

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 e^–S_I, 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 ^ze^–wt¢. (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 U_t (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 (^¶U_t /^¶z)_{z = 1} is for n = 0 (the sink itself), which is ze^z e^–wt.

For r = 2 the sink has two lines. There is then the n = 0 contribution, z² e^z e^–2wt, and the 1–vertex diagram shown in Fig. 3. Including the combinatorial factor of 2, its contribution is

Using these results, we readily find that

For n(0) = n₀, 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 án^rñ_f.

7 Malthus-Verhulst Process

We now return to the Malthus-Verhulst process; the three vertices associated with –S_I [see Eq. (54)] are shown in Fig. 4. The series for U (z) begins with the "sink" term, z^r e^z e^–rwt, and then includes diagrams with n > 1 vertices. (Note that the n-th order contributions carry factors of (–n)^s

^n–s, with s variable.) Based on our diagrammatic analysis of the birth-and-death process, we can formulate the following simple rules for evaluating the n–th 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 t₁,...,t_n to the vertices from left (nearest the sink) to right.

3) For each external line include a factor of ze–^w, 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 e^{–w(t_i – t_j)}.

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 t₁,...,t_n, with

(Note that by fixing the time ordering we effectively include the factor 1/n! that comes from expanding e^–S_I.)

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 án^r ñ_f:

We have seen that each line entering a node with time variable carries a factor of e^–w, and each line exiting the node carries a factor of e^wt. 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 n₀e^–wt may be written as

with A = n(n₀ –1) and B = n[(n₀ –1)(n₀ –3)n– ] . Evidently, our analysis generates an expansion in powers of t, whose convergence would seem to require that nn₀ 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 – e^–wt 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/w² To proceed, we form a Padé approximant [16,17] to the power series in y:

Equating coefficients of y and of y², 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) ñ = n₀ e^–wtf(t) (with f represented by the Padé approximant), and with the simple exponential decay, n₀ e^–wt, for n₀ = 3 and n₀ = 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 n₀ = 10 we have nn₀ = 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 n_r(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 U_t is a function of the variables z_r and z_r defined at each site; we use {z} and {z} to denote these sets of variables. To write U_t for this process we add the diffusive contribution to the Malthus-Verhulst part, Eq. (51), to find

where Dy_r = å_e (y_r+e – y_r) is the discrete Laplacian, and the functional integrals are now understood to include the funcions y_r 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 U_t 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 y_r(t)® y(x, t), and Dy_r ® Ѳ y.

iii) We discard the term µ

² y², 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 d_c = 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 µ y² . 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, – y² , 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 w_c = 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 µ

ⁿ, 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,b–diagram has a weight of 2^n–b–1. Let W(n,b) be the sum the weights of all n,b–diagrams. 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,b–diagram, we can generate a set of distinct n+1–vertex 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 n–vertex labelled diagrams are mutually disjoint; (3) applied to the complete set of n–vertex labelled diagrams, the procedure generates the complete set of n+1–vertex 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,b–diagram has b+1 vertices of degree zero, and n–2b–1 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 xⁿ y^b/n! and sum over n and b. Letting m = n–1 and rearranging, one finds

Now let b¢ = b–1 in the terms multiplying W(m,b–1). After a simple rearrangement (and dropping the primes) we have

We require

For y = 1, Eq. (134) immediately yields G(x) = x/(1–x), 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

Appendix

), 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 ln_t 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 n–2b–1 nonterminal vertices that are not bifurcations, which carry factors b_j. Each b_k carries a factor of e^–wt (z–1) for a total of n+1 such. Next consider the factors e^wt_i for each vertex i. If i is terminal, there are two such factors (from b), and a factor e^–wt_i due to the line entering the vertex. If i < n is neither terminal nor a bifurcation, there is a factor e^wt_i from b_i, while the exponential factors associated with the lines entering and leaving the vertex cancel. If i < n is a bifurcation, there is no factor b_i, but there is again a net factor of e^wt_i as two lines exit, while only one enters vertex i. Consider vertex n. If it is a bifurcation, then the net factor is e^wt_n (two lines exiting, factor e^–wt_n in c_n). The same holds if vertex n is not a bifurcation (single line exiting, factors in b_i and c_i cancel). In summary, there is a factor e^wt_i associated with each vertex. Finally, all vertices carry a factor of 2, except for the terminal ones, leading to an overall factor of 2^n–b–1. Combining all of these observations, the contribution due to a given labelled n,b–diagram is

Address to correspondence

Ronald Dickman, Ronaldo Vidigal

Departamento de Física, ICEx, Universidade Federal de Minas Gerais,

30123-970, Belo Horizonte, MG, Brasil

Publication Dates

• Publication in this collection
  23 Apr 2003
• Date of issue
  Mar 2003

History

• Received
  17 Aug 2002