TRAVELLING WAVES FOR A CERTAIN FIRST-ORDER COUPLED PDE SYSTEM

. This paper concentrates on a particular (cid:12)rst-order coupled PDE system. It provides both a detailed treatment of the existence and uniqueness of monotone travelling waves to various equilibria, by di(cid:11)erential-equation theory and by probability theory and a treatment of the corresponding hyperbolic initial-value problem, by analytic methods. The initial-value problem is studied using characteristics to show existence and uniqueness of a bounded solution for bounded initial data (subject to certain smoothness conditions). The concept of weak solutions to partial di(cid:11)erential equations is used to rigorously examine bounded initial data with jump discontinuities. For the travelling wave problem the di(cid:11)erential-equation treatment makes use of a shooting argument and explicit calculations of the eigenvectors of stability matrices. The probabilistic treatment is careful in its proofs of martingale (as opposed to merely local-martingale) properties. A modern change-of-measure technique is used to obtain the best lower bound on the speed of the monotone travelling wave | with Heaviside initial conditions the solution converges to an approximate travelling wave of that speed (the solution tends to one ahead of the wave-front and to zero behind it). Waves to di(cid:11)erent equilibria are shown to be related by Doob h -transforms. Large-deviation theory provides heuristic links between alternative descriptions of minimum wave speeds, rigorous algebraic proofs of which are provided.


Introduction
Connections between parabolic and elliptic equations and diffusion processes are well-known. In this paper we explore the connections between hyperbolic partial differential equations and stochastic processes by treating the following system. Let r 1 , r 2 , q 1 , q 2 be positive (since they correspond to rates of probabilistic processes) constants and b 1 , b 2 unrestricted real constants, fixed throughout. Let θ be a positive rate parameter. We consider an equation system related to the generalized FKPP system discussed in Champneys, Harris, Toland, Warren and Williams [2]. Here the system of interest is where u is a vector-valued function from [0, ∞) × R to R 2 , u = (u 1 , u 2 ) ∈ R 2 and u 2 = (u 2 1 , u 2 2 ), and where (We use ':=' to mean 'is defined to equal'.) Dunbar [8] considered a similar system but with q 1 = q 2 , r 1 = r 2 and the nonlinearity was 1 4 (u 1 + u 2 ) 2 1 1 rather than u 2 . Holmes [10] compared reaction-diffusion systems to reaction-telegraph models (again with symmetry in the nonlinearity in both terms) for animal movement; if equation (1) is rewritten as a single second-order PDE it can be seen that our system is a reaction-telegraph model, though we find it more convenient to work with the pair of first-order equations. Hadeler [9] discussed more general non-linearities but still retained the condition q 1 = q 2 . We contrast the probabilistic interpretation of these models and ours in section 1.4. A travelling-wave solution of equation (1) The 'source point' S = (0, 0) and the 'target point' T = (1, 1) are clearly equilibria of equation (2). If r 1 r 2 ≥ 4θ 2 q 1 q 2 , there will also be equilibria at the two points where ρ i := q i /r i and ∆ := 1 4 − θ 2 ρ 1 ρ 2 . For θ ∈ (0, ∞) we study the existence of monotone travelling waves from S to T . Waves from S to the other equilibria can be obtained through a transformation detailed in section 1.3. Monotone travelling waves from S to T have a direct probabilistic interpretation which we will explore in this paper, but whether such an intepretation of other travelling waves is possible is unclear.

Stability of equilibria
Suppose that (B + cI) is invertible and write equation (2) in the form Then F is a quadratic polynomial. Let E be an equilibrium point of (2), thus F (E) = 0. Then write w(x) − E = v(x) and expand equation (3) to first order in v. This yields Thus K c,θ (T ) (we write K c,θ to emphasize the dependence on c and θ) satisfies (B + cI)K c,θ (T ) + R + θQ = 0 (4) and K c,θ (S) satisfies (B + cI)K c,θ (S) − R + θQ = 0.
The stability properties of these matrices are investigated in section 3.

The main ODE theorem
DEFINITION. An eigenvalue λ of a real 2 × 2 matrix M will be called stable [respectively, unstable] monotone if This definition links nicely with the Perron-Frobenius Theorem (see Seneta [15]). The theorem implies that a square matrix M with all off-diagonal entries strictly positive has a special eigenvalue Λ P F (M ) with an associated eigenvector with all entries positive and such that every other eigenvalue of M has real part less than Λ P F (M ). Moreover, any eigenvector with all entries positive must be a multiple of the Perron-Frobenius eigenvector.
As will be proven in section 3 (Lemma 3.1), for any fixed θ > 0, there exists a critical value c(θ) in the interval min then K c,θ (T ) has at least one stable monotone eigenvalue, and if c < c(θ), c = min(−b 1 , −b 2 ), it has no such eigenvalues. Since a necessary condition for the existence of a monotone travelling wave of (2) which converges to T as x → ∞ is the existence of a stable monotone eigenvalue of K c,θ (T ) this provides a lower bound on possible values of c for which a monotone travelling wave can exist. It is shown later that this condition on c is also sufficient for the existence of a monotone connection from S to T , see section 4.5.
Large-deviation theory gives probabilistic heuristics for this critical value, see section 5.2. Thus we can state our main result as follows. This theorem is proven in section 4 using a shooting argument. Behaviour at c(θ) itself depends on the relative values of parameters as follows (these special cases are detailed in section 4.6):

A Doob h-transform
As in Champneys et al. [2] there is a Doob h-transform that maps monotone waves from S to one of the equilibria to monotone waves from S to another of the equilibria. This fact is summarised by the following Lemma.

Lemma 1.2
Suppose that θ is fixed at a value where ∆ ≥ 0, so that E + and E − exist. If E = (α 1 , α 2 ) denotes either E + or E − , then the substitutioñ q i := q i α j /α i (j = i),r i := r i α i ,ũ i := u i /α i ,w i := w i /α i (5) transforms (1) and (2) into their ∼ versions, monotone waves from S to E for the original problem corresponding exactly to monotone waves from S to T for the ∼ problem. The possibility that E + = E − and ∆ = 0 is not excluded.
Much of our work therefore automatically transfers to the case when T is replaced by E + or E − , though the critical values c ± (θ) of θ corresponding to waves from S to E ± will be different.

Key probability theorems
Our aim in section 5 will be to prove the following theorems which connect a probability model, defined below, with the system (1). Note that we often switch between two equivalent notations in the probability part of this paper to reduce the use of subscripts. Thus we will sometimes write, for example, b(y) for b y , w(x, y) for w y (x) and u(t, x, y) for u y (t, x) (for y = 1, 2).
Let I := {1, 2} and consider the following two-type branching system of particles. At time t ≥ 0, there are N (t) particles, the k-th particle -in order of birth -having position X k (t) in R and type Y k (t) in I. The state of the system at time t is therefore Particles, once born, behave independently of one another. Each particle lives forever. The type of a particle (once born) is an autonomous Markov chain on I with Q-matrix θQ. While a particle is of type y ∈ I, it moves with constant velocity b(y), and it gives birth -to one child each time, at its own current position and of its own current type -in a Poisson process of rate r(y). So, r(y) is the breeding rate of type y.
The branching system in Dunbar [8] had q 1 = q 2 , r(1) = r (2) and particles, rather than giving birth to one particle of the same type and living on themselves, die and give birth to two particles of independent random types, with each having equal probability of being of either type. The equations studied by Hadeler [9] correspond to the new pair of particles having type correlated to that of their parents, for any correlation except ±1. Our model permits q 1 = q 2 , so that the particle-type Markov Chain can have any equilbrium distribution on I For our model it makes no difference whether you consider that one new particle has been born and the old one lives on too, or that two new particles of the same type are born. However, when the type of the new particles is random, this distinction is crucial. Our model is thus distinct from those previously considered.
Write P x,y (with associated expectation E x,y ) for the law of this process when it starts from one particle of type Y 1 (0) = y at position X 1 (0) = x. By martingale [respectively, local martingale, supermartingale, ...] we mean a process which is for every P x,y a martingale [respectively, ...] relative to the natural filtration F t (P x,y -augmented, to be precise) of the process at (6).
The state-space for this process is Define L(t) := inf k≤N (t) X k (t). This is the position of the left-most particle. The asymptotic speed of the left-most particle is lim t→∞ t −1 L(t).

Theorem 1.3
As t → ∞, the following holds almost surely (a.s.) If u satisfies the coupled system (1), for t ≥ 0 and x ∈ R and if , and u is an approximate travelling wave of speed c(θ) in the sense that This theorem is proved, along with the following one, in sections 5.4 and 5.5. The theorem allows us to relate the speed of the spread of the particles in the probabilistic model with the wave speed of travelling waves. Specifically, we are claiming that the left-most particle travels, in the limit, at the speed −c(θ), where c(θ) is the critical speed above which a unique monotone travelling wave exists and below which no such wave exists.
An analytic proof that the weak solution of (1) is between 0 and 1 for the Heaviside initial data is included in section 2.4, as well as a proof that for continuous initial data between 0 and 1 the solution remains between 0 and 1 (see Lemma 2.2 and the subsequent remarks). The work of section 2 is not necessary for the probabilistic approach to the problem but adds insight from the viewpoint of classical analysis, and vice-versa. In section 5 we show that any (smooth) solution to the coupled system (1) that is between 0 and 1 has a McKean representation. We use this representation to motivate the (probabilistic) construction of a solution for the initial-value problem with Heaviside initial data. This constructed solution is then directly verified to satisfy the appropriate equations and does remain between 0 and 1.

This implies that
so that (8) follows.
This result gives us the weaponry to prove Theorem 1.3 using martingale techniques. Note that Theorem 1.3 is proved using Theorem 1.4, and that in section 5.4 we prove Theorem 1.4 via the Theorems and Lemmas in the preceding sections of section 5 -hence various references to Theorems of section 5 in the statement of Theorem 1.4 are not circular.

Summary chart of results
This  (1) becomes an ODE and an algebraic equation. The only candidate for a travelling wave from S to T is the corresponding segment of the solution of the algebraic equation (which is a parabola) -full details of this special case are in section 4.6.

Existence and uniqueness results for the PDE initial-value problem
To study the question of global existence and uniqueness of solutions to the initial-value problem for the system (1), it is convenient to change to moving coordinates (moving at a speed of 1 2 (b 1 + b 2 )) and then re-scale space so that the coefficients of u x are 1 and −1. This is possible unless b 1 = b 2 -we deal with this case in section 2.6. We use subscript notation to represent derivatives and relabel so that the system becomes: The functions f and g are introduced to simplify notation.
We are particularly interested in the Cauchy problem for Heaviside initial data: which we study in section 2.4, but first we study the Cauchy problem with C 1 -initial data.
The characteristics are x+t = constant, on which ∂u ∂t = f (u, v) and x−t = constant, on which ∂v ∂t = g(u, v). Integrating along the characteristics from t = 0 to a point (T, X) gives: We work in the Banach space C 1 of 2-vector-valued functions w for which w and w x are continuous and bounded for all x, in which we choose as norm |||w||| = max(||w||, ||w x ||), where ||w|| := sup x∈R |w(x)|, the usual L ∞ -norm. These equations can be used with the Contraction Mapping Principle to prove existence and uniqueness of classical solutions for a short time for C 1 -initial data. See Courant and Hilbert [5, pages 461-471] for details.

Proofs for smooth initial data
To go from short-time existence and uniqueness of solutions to global existence and uniqueness we prove the following lemmas which give bounds on the solutions for all time (for a certain class of initial data). These bounds then allow iteration of the Contraction Mapping argument -hence local existence and uniqueness become global. This iteration is done by using the solution obtained from local existence and uniqueness, up until a small, fixed time, τ , then taking the value of this solution at time τ as initial data, and repeating the argument. The bounds obtained in the following two lemmas allow repeated use of the same τ at each step, rather than having to take a sequence of τ n (whose sum may converge to a finite blow-up time), thus we obtain existence and uniqueness for all time. To set our notation, note that we will write u(0, (9), (10) satisfy 0 < u, v < K for all time.
Proof. Suppose, for a contradiction, that there exists a point (T, X) where (u, v) is outside the square, (0, K) 2 . The value of (u, v) at this point only depends on the initial data in the interval [X − T, X + T ]this is the domain of dependence (see Courant and Hilbert [5,). This data determines the solution (u, v) throughout the closed triangle, which we shall denote by Ω, of (t, x)-space whose corners are (T, X), (0, X − T ) and (0, X + T ).
Since Ω is compact and u and v are continuous, for each of the possible violations (that is, violations of the four inequalities u < K, u > 0, v < K and v > 0) we can find a first time it occurs in Ω. For example, if u ≥ K at some point in Ω then there exists a point (t 0 , x 0 ) ∈ Ω such that u(t 0 , x 0 ) = K and, for all (t, x) ∈ Ω with t < t 0 , u(t, x) < K. Similarly we can find a time (and corresponding spatial position) where the first violations of u > 0, v < K and v > 0 occur in Ω (if such violations do occur). We can then study the first violation that happens, by taking the one corresponding to the minimum of these 4 times (taking the time to be T + 1 if it does not occur in Ω). This is well-defined because we know there is at least 1 violation and at most 4.
So, consider each possible first violation in turn. Firstly, that there exists a point (t 0 , x 0 ) ∈ Ω such that u(t 0 , x 0 ) = K and, for all (t, x) ∈ Ω with t < t 0 , u(t, x), v(t, x) ∈ (0, K) 2 . Then consider u restricted to the characteristic x = x 0 + t 0 − t through (t 0 , x 0 ) (which lies entirely in Ω). From (9), with u = u(t, x 0 + t 0 − t), Similarly, we can show that v(t, x) < K for all (t, x) ∈ Ω, and we use a similar argument to show u > 0 and v > 0: Hence no violation occurs in Ω which is a contradiction. This completes the proof.
Lemma 2.1 allows us to prove the following: Proof. Consider a sequence (u n , v n )(x) of initial data satisfying the conditions of Lemma 2.1, which converges in This lemma completes the proof of global existence and uniqueness of bounded solutions for smooth initial conditions between 0 and any constant 0 < K ≤ 1. This upper limit on K is the best possible, since u 2 − u > 0 for u > 1, so the solution tends to grow. Indeed, it is clear that for u 0 (x) = v 0 (x) identically equal to 1 + , for any > 0, the solution blows up in finite time.
However it is possible to extend the above lemmas to deal with initial data below 0, since the u 2 − u nonlinearity will tend to push the solution up towards 0.
Proof. We can again look for the first violation, this time of the four restrictions u > K, u < 0, v > K, v < 0, and look at a characteristic going through a point at which the first violation occurs.
For the case u = K, note that so that u does not hit K.
For u = 0, note that, looking at the characteristic sufficiently close to the violation point (so that u > −1) so that u does not hit 0 either.
Similarly for v. Thus all the contraction mapping arguments extend to all time in the same way that the bounds on data between 0 and K (for 0 < K ≤ 1) extended existence and uniqueness in that case.
In the same way that we passed from Lemma 2.1 to Lemma 2.2 we can relax the strict inequalities in the hypotheses of Lemma 2.3 to obtain the following result.
Combining Lemma 2.2 and Lemma 2.4 prepares the ground for the following result.
Proof. If K 1 = 0 then this result is simply a restatement of Lemma 2.2, and Lemma 2.4 deals with the case K 2 = 0. So, we may assume that For a violation such as, say, u = K 1 , note that, considering a section of the characteristic sufficiently close to the violation point for u < 0 so that u does not hit K 1 .
For u = K 2 , note that, considering the characteristic sufficiently close to the violation point (so that u > 0) Similarly for v, hence we are done.
Finally, again using the inequality relaxation, and noting that if the initial data is identically zero then the solution is identically zero for all time, we have: No such result will be true for a pair of constants K 1 and K 2 both strictly on the same side of zerothe solution to the initial value problem with u 0 = v 0 identically equal to some constant K, such that −∞ < K < 1, tends monotonically to zero (see section 2.6). This fact is put together with comparison arguments in section 2.2 to give a much stronger result than Lemma 2.6.
Given bounded initial data with an upper bound no greater than 1, we can read off the appropriate values of K 1 and K 2 by defining: Local existence of solutions with piecewise smooth initial data (i.e. continuous but with a finite number of jump discontinuities in the x-derivative) follows by approximating (in the supremum norm) piecewise smooth data by C 1 data -the limiting solutions are classical except on characteristics x ± t = x 0 propagating from points x 0 of discontinuity of ∂u0 ∂x , ∂v0 ∂x .

Comparison results
Consider two solutions (u, v) and (ũ,ṽ) of the equations (9),(10), with u(0, On the characteristics of the form x + t = constant, we know that and similar equations for v andṽ along the other characteristics.
Studying the difference between the solutions along a characteristic x + t = constant, we see that it obeys the equation Provided that the initial data (u 0 , v 0 ) and (ũ 0 ,ṽ 0 ) satisfies the conditions of Lemma 2.6, then (u +ũ) will be bounded below (for all time), enabling us to write, along the characteristic up until a putative equality of u andũ, for some constant C. There is a similar equation for the difference of v andṽ.
Thus, if u 0 >ũ 0 and v 0 >ṽ 0 , then, for all time, u >ũ and v >ṽ. By taking sequences of initial data we obtain the following result.
Since the solution to the initial value problem with u 0 = v 0 identically equal to some constant K, such that −∞ < K < 1, tends monotonically to zero (see section 2.6), then any solution bounded between −∞ < K 1 ≤ 0 and 0 ≤ K 2 < 1 will tend to zero. Thus Lemma 2.6 can be modified to the following result.

Weak solutions and the Rankine-Hugoniot jump conditions
To deal with discontinuous initial data it is necessary to utilize the concept of weak solution.

10
DEFINITION. For u, v bounded and measurable, we say that U := (u, v) is a weak solution of equation (1) if, for all test functions Φ, For a piecewise classical solution to be a weak solution its curves of discontinuity in the xt-plane must satisfy certain conditions, known as the jump or Rankine-Hugoniot conditions. For a semi-linear system such as (1) it is well known that these conditions imply that a component of the weak solution can only have discontinuities across characteristics corresponding to that component (shock paths can only be characteristic curves).
Moving to the weak setting requires care since jump conditions are not necessarily sufficient to guarantee uniqueness of solutions. Lemma 2.9 proves existence and uniqueness of weak solutions for a certain class of bounded, measurable initial data.

Heaviside initial data
We now construct explicitly a piecewise classical solution for Heaviside initial data that satisfies the Rankine-Hugoniot jump conditions (hence it is a weak solution and will be shown to be unique by Lemma 2.9).
For the Heaviside initial data the jump conditions reduce to two requirements -that the discontinuity in u propagates along the characteristic for u that goes through zero, i.e. u jumps across x = −t and is continuous elsewhere, and that the discontinuity in v propagates along the characteristic for v that goes through zero, i.e. v jumps across x = t and is continuous elsewhere. We have defined the Heaviside function so as to be left-continuous and construct so that our solution inherits this property -this matches the continuity which the probabilistic approach implies. The probabilistic interpretation of this solution is discussed in section 5.5.
Clearly the solution for the Heaviside data is identically zero for x ≤ −t and identically one for x > +t.
We are constructing a left-continuous solution -therefore on x = +t, u = 1 while v remains to be calculated.
We integrate along the x = +t characteristic to find v there. Since u = 1, and v(0, 0) = 0. Hence, when and v increases from 0 to min θq2 r2 , 1 as t goes from 0 to ∞. When and v increases from 0 to 1 as t goes from 0 to ∞.
We can also investigate the discontinuity in u along x = −t similarly. v is continuous across this characteristic and zero on it, so is zero on x = −t+. Thus, integrating along the inside edge of the characteristic, and u(0, 0+) = 1. Hence, and u decreases from 1 to 0 as t goes from 0 to ∞.
We now know the values of both u and v on the inside edge of the wedge |x| < t, between the discontinuities -and this data is continuous and between 0 and 1. Thus there is a unique solution to the problem with these as initial/boundary values and this solution is between 0 and 1. It can be found by following characteristics x + t = constant from the x = +t characteristic. This solution is piecewise classical with the only discontinuities being across the characteristics as required and so is indeed a weak solution.
Results of computational work on this initial-value problem will be reported elsewhere.

Discontinuous initial data
We follow the method of Beale [1] in his work on the Broadwell model. Let φ be a test function with For non-smooth initial data (u 0 , v 0 ) we consider initial data formed by mollifying (u 0 , v 0 ). We use the convolutions u In the following argument we consider differences between solutions for different mollifications of the same measurable initial data. We work with initial data that is bounded between constants K 1 and K 2 , For the argument to work these differences between solutions should be in L 1 , which we show follows if the difference between two different mollifications of the initial data is in L 1 . This is clearly true for the Heaviside data, and for initial data that is itself in L 1 . Then, for example, it is true for bounded initial data that only differs from the Heaviside function by an L 1 -function.
With this in mind define a step function s to be a function of the form for some constants k 1 and k 2 such that −∞ < k 1 , k 2 ≤ 1. Then if both components of the initial data differ from step functions (with possibly different constants) by L 1 -functions, the difference between mollifications will be L 1 .
. Then the u (k) and v (k) are between K 1 and K 2 for all t, x and, for 1 ≤ p < ∞, for τ > 0 arbitrary, . This weak solution is also bounded by K 1 and K 2 and is unique.
Multiplying equation (12) by f n (U i ), yields: Integrating in x over a finite interval [x 1 , Letting n → ∞ and taking limits using the bounded convergence theorem implies that where in the latter step we are using the fact that the initial data is only an L 1 -function from a step function to say that Choosing T sufficiently small so that CT < 1 and rearranging gives, where the constant on the right-hand side (which we will now denote by γ) does not depend on x 1 and x 2 (γ is finite due to boundedness of solutions for smooth, bounded initial data). Thus, for 0 ≤ t < T , and hence each U i (t, ·) ∈ L 1 ((0, T ) × R). Now it is possible to return to equation (13) and improve our estimate. Firstly, rewrite it as: Since each U i (t, ·) ∈ L 1 ((0, T ) × R) there exist sequences such that, letting x 1 → −∞ and x 2 → ∞ along these sequences, keeping t fixed, the latter two terms of equation (14) tend to zero. The fact that each Hence: Thus, from Gronwall's inequality, . In fact, since the u (k) and v (k) are all bounded between K 1 and K 2 , the convergence takes place in C 0, τ; is a weak solution of the equation and that (u, v) is again bounded by K 1 and K 2 for all time.
We now verify the uniqueness property. For these solutions we can show, for u that, for arbitrary and similarly for v. Therefore for almost all (t, x). Hence if there are two solutions with the same initial data and y(t) is the L ∞ -norm of the difference at time t, we obtain (using the boundedness between K 1 and K 2 ) an estimate and this implies y is identically zero.
The restriction to initial data that is an L 1 function different from a step function can be removed using domains of dependence and a truncation argument as follows.
The value of u on the line (T, X) to (T, X + 1) depends only on the initial data in the interval [X − T, X + T + 1]. Thus it should agree with the solution whose initial data is identically zero outside this interval, and matches on the interval, which we shall describe as the truncated initial data. However the truncated initial data satisfies the conditions of Lemma 2.9 and so the truncated initial data has a unique corresponding solution. Thus we can define the solution for more general initial data to be that constructed by piecing together solutions for the truncated initial data, Lemma 2.9 guarantees that this will be well-defined. 14 2.6 The case when b 1 = b 2 When b 1 = b 2 we change to moving coordinates at speed b 1 . Relabelling as before we obtain the ODE: For this pair of equations the square [0, 1] 2 is positively invariant -since f (0, v) > 0 and f (1, v) < 0 for 0 < v < 1 and similarly g(u, 0) > 0 and g(u, 1) < 0 for 0 < u < 1. S and T are clearly equilibria. At all other boundary points the flow is into the unit square in forwards time. In fact, any square of the form For the Heaviside initial data the solution is clear -the Heaviside function simply propagates at speed −b 1 . This is also clear from the probabilistic method (see section 5.5).
The pair of equations (15) is also relevant because, if the initial data in the PDE system (equation (1)) The constant initial data case is particularly relevant in light of the comparison arguments in section 2.2. Solutions with initial data bounded below and above by constants K 1 and K 2 can be bounded below and above by the solutions for constant initial data K 1 and K 2 , but these solutions clearly go to zero for −∞ < K 1 ≤ 0 and 0 ≤ K 2 < 1 (since the square [K 1 , K 2 ] 2 is positively invariant, contains only one equilibrium, which is S, and a periodic solution within the square is not possible due to the direction of the vector field).

Useful algebraic results
Throughout this section we shall assume, without loss of generality, that b 1 ≥ b 2 . Proof. In fact, we can say rather more than this about K c,θ (T ). Explicitly, if (B + cI) is invertible, Thus, its two eigenvalues, λ + and λ − , are given by: These eigenvalues correspond to eigenvectors of the form Bifurcation diagrams of the eigenvalues as functions of c are given in Figures 1 and 2. For c > −b 2 , note that λ − is a stable monotone eigenvalue and λ + is not, since v + < 0 (as λ + > θq1−r1 b1+c in this case). For c < −b 1 , note that neither eigenvalue is stable monotone as v − < 0 and λ + > 0. Thus, if there is a critical value c(θ) it lies in the interval claimed. For b 1 = b 2 = b, say, this is enough to prove the lemma, c(θ) = −b, independent of θ.
So, to complete the proof of the lemma it is sufficient to consider the case b 1 > b 2 , where it is necessary to examine the eigenvalues when −b 1 Examining the term under the square root in equation (16) (which we will denote by h(c), for fixed θ), for −b 1 < c < −b 2 , note that it can have one, two or no zeroes. If there are two zeroes let us denote them by c 1 and c 2 , with −b 1 < c 1 < c 2 < −b 2 . When θq 2 > r 2 and θq 1 = r 1 , there are two zeroes of h and for c 1 < c < c 2 , h(c) < 0. So in this case the eigenvalues of K c,θ (T ) are complex. When −b 1 < c < c 1 and when c 2 < c < −b 2 , h(c) > 0, so the eigenvalues are real. When c 2 < c < −b 2 both eigenvalues are stable monotone. When −b 1 < c < c 1 the eigenvalues are not stable monotone -θq 1 < r 1 implies that the eigenvectors have components of opposite signs (as shown in Figure 1 (a)) while θq 1 > r 1 implies that the eigenvalues are positive (as shown in Figure 1 (b)). Thus, in this case, c(θ) = c 2 .
When θq 2 > r 2 and θq 1 = r 1 there is a single zero of h, at c 1 say (as shown in Figure 1  Hence the lemma is proven. When c = c(θ) there is a repeated stable monotone eigenvalue, unless c(θ) = −b i , in which case the point is moot since the stability matrix does not exist, but the analysis can be completed by direct methods (see section 4.6).
For the probabilistic method we need to look at one of these special cases in more detail -the case Recall that K c,θ (T ) is defined by equation (4), thus an eigenvector v (with corresponding eigenvalue λ) of K c,θ (T ) will satisfy the equation For (B + cI) invertible, this relation is also true in the opposite direction -a non-trivial vector v that satisfies equation (18) is an eigenvector of K c,θ (T ) with eigenvalue λ. However, for (B + cI) singular, equation (18) can still have non-trivial solutions, in the case of interest it has one solution, Thus this solution (which is what we are really interested in, analysis of K c,θ (T ) is a short-cut, and simplifies discussion by giving us a way of referring to eigenvalues and eigenvectors as those of K c,θ (T )) can be thought of as a generalized stable monotone eigenvalue of K c,θ (T ) and we have proved that, for c > c(θ) the matrix K c,θ (T ) has at least one (possibly generalized) stable monotone eigenvalue.
The value of c(θ) can also be derived from the probabilistic model of the system: the theory of large deviations gives a formula for c(θ) which effectively summarises all these cases; see section 5.2. That formula, which is shown in section 5.2 to be entirely equivalent to the preceding characterization, makes it easy to observe that c(θ) is decreasing as θ increases and that This limit will be discussed in section 5.2, for now note that the lower bound on c(θ) given by Lemma 3.1 can therefore be tightened to c * .
It is also possible to arrive at this limit by further manipulating h(c). Rearranging the equation h(c) = 0 to express it as a quadratic in c (where the coefficients are quadratic in the other parameters, including θ) enables us to consider the limit as θ → ∞ by picking out the terms of highest order (i.e. order two in θ) only. Simplifying this expression yields The following result is required in our proof of L 1 convergence of Z λ in Theorem 5.4.

Lemma 3.2 (i) Suppose that c > c(θ).
Let λ s (c) be the stable monotone eigenvalue of K c,θ (T ) (the one nearer to 0 if there are two). From the definition of K c,θ and λ s (c) it is easily seen that −λ s (c)c is the Perron-Frobenius eigenvalue of (λ s (c)B + θQ + R). For µ < λ s (c), with µ sufficiently close to λ s (c), (iii) When c(θ) < −b 2 , K c(θ),θ (T ) has a double eigenvalue, which we will denote by λ 0 . This eigenvalue is geometrically simple, i.e. it has only one normalised eigenvector even though it has algebraic multiplicity two, and is stable monotone.
(ii) The first part follows from the fact that λ s (c) −1 Λ P F (λ s (c)B + θQ + R) = c. For the second part consider again the explicit form of the Perron-Frobenius eigenvalue as a function of µ: It is clear, by continuity (considering Λ P F as a function of µ and noting that the term under the square root is strictly positive), that We claim that Λ P F (µB + µc(θ)I + θQ + R) must attain a local minimum at µ = λ 0 . Note that, for fixed c > c(θ), sufficiently close to c(θ), there are two values of µ such that Λ P F (µB + µcI + θQ + R) = 0, these two values of µ are the two stable monotone eigenvalues of K c,θ (T ). We now use a convexity result due to Cohen [4].
To be precise, Cohen's result states that the Perron-Frobenius eigenvalue of a matrix M 1 + M 2 is a convex function of M 2 , where M 1 has positive elements off the main diagonal and M 2 is a diagonal matrix (possibly 0). Thus, for µ between the two eigenvalues Λ P F (µB + µcI + θQ + R) < 0 (strict inequality since there can be at most two zeroes -each is an eigenvalue of K c,θ (T ), a 2 × 2 matrix). Thus, for c = c(θ), Λ P F (µB + µc(θ)I + θQ + R) > 0 except at µ = λ 0 where it is zero.
Thus the derivative on the right hand side of equation (19) is zero at µ = λ 0 and hence the proof of (ii) is complete.
(iii) As c decreases through c(θ) the two stable monotone eigenvalues of K c,θ coalesce and, at least for c sufficiently close to c(θ), become a complex conjugate pair with negative real part. From equation (16), b1+c(θ) + θq2−r2 b2+c(θ) , and from equation (17), this corresponds to an eigenvector of the form We can verify directly that λ 0 and v 0 as given above have the correct signs (negative and positive respectively) by using inequalities arising from the definition of c(θ).
For the probabilistic proof of uniqueness, modulo translation, of monotone travelling waves from S to T , the next lemma is important.

Lemma 3.3 Suppose that c > c(θ)
and that K c,θ (T ) has two stable monotone eigenvalues. Let β be the stable monotone eigenvalue further from 0. Then, for α > β with α sufficiently close to β, the only non-negative 2-vector g such that 0 ≤ α(B + cI) + θQ + R g is the zero vector: g = 0.
Proof. Convexity of Λ P F (µ(B + cI) + θQ + R) as a function of µ given by Cohen's Theorem [4], for β < µ < λ, yields Λ P F (µ(B + cI) + θQ + R) ≤ 0. However, there can only be two values (µ 1 and µ 2 , say) of µ at which this function is 0 since each µ i is thus an eigenvalue of a 2 × 2 matrix. These two values are therefore β and λ and the inequality is strict.

Plan of attack
To prove Theorem 1.1 we use shooting arguments. First note that the path of any solution w of (2) which satisfies must lie entirely inside the open unit square. Whether it must also lie between the nullclines w 1 = 0 and w 2 = 0 depends on the values of various parameters, as explained below. Depending on the geometric configuration of the nullclines, we introduce shooting boxes. We observe that these regions have four key properties upon which our proof will be based: • Non-constant solution curves do not intersect the boundary of the region tangentially and hence any non-constant solution curve which intersects the boundary crosses it transversally; • Exit-times of solutions in the regions are continuous functions of initial conditions; • No non-constant solution curve passes through S or T ; • Other than S and T there are no equilibria in the closure of these regions.
There are two types of region that we use -inside regions when c > max(−b 1 , −b 2 ) and outside regions when min . The terminology is based on the geometry of the phase plane as will be seen shortly.
Also note that the relative position of T and the maximum of the parabola P i is determined by the sign of θq i − r i -if θq i ≥ r i then the segment of P i from S to T is monotone; if θq i < r i then this segment has a turning point before reaching T .
Regions of interest. Now letΣ denote the rectangle in the w-plane with two sides on the axes intersecting at 0, a side through E 1 parallel to {w 1 = 0} and one through E 3 parallel to {w 2 = 0}. Let Σ denote the open convex subset ofΣ whose boundary comprises four straight-line segments from ∂Σ, a parabolic segment from P 1 joining E 1 to E 2 and a parabolic segment from P 2 joining E 2 to E 3 . Thus the boundary of Σ always has four straight-line segments: in addition it has two parabolic components when E + , E − and T are distinct, one parabolic component when two of E + , E − and T coincide and no parabolic component when θ > (4ρ 1 ρ 2 ) − 1 2 . Also Σ consists of the union of three sets: ω 2 =Σ ∩ Ω 1 ∩ Ω 2 , a relatively closed component ω 3 whose boundary intersects {w 1 = 0} away from the origin and a relatively closed component ω 1 whose boundary intersects {w 2 = 0} away from the origin. (See Figure 3.) Note that and hence, if w = (w 1 , w 2 ) satisfies (2) then, A monotone curve connecting S to E 1 must approach E 1 from ω 1 ∪ ω 2 . Thus, by (20), if (b 1 + c) < 0 there cannot be a monotone connection from S to E 1 . Similarly, a monotone curve connecting S to E 3 must approach E 3 from ω 2 ∪ ω 3 . Thus, by (20), if (b 2 + c) < 0 there cannot be a monotone connection from S to E 3 .
A monotone connection to E 2 must eventually approach through ω 2 , so a necessary condition for existence is (b 1 + c) > 0 and (b 2 + c) > 0.
When considering monotone connections from S to T it is sufficient to restrict attention to the intersections of each ω i with the unit square -a monotone connection to T must lie entirely within the unit square and must eventually approach T through some ω i . It is possible to work with just ω 2 and ω 3 -since ω 1 is mapped to ω 3 by swapping b 1 and b 2 , r 1 and r 2 and q 1 and q 2 . (It is not possible for a monotone connection from S to T to lie partly within ω 1 and partly within ω 3 since both w 1 and w 2 have signs opposite in ω 1 to their signs in ω 3 .) The outside region is ω 3 intersected with the unit square; the inside region is a subset of ω 2 inside the unit square. Examination of these two regions allows us to cover all possible monotone connections to T -connections must eventually approach from one of the two.
Also note that the nullclines represent points at which solution curves are vertical or horizontal (where the nullclines cross, there are fixed points). Thus, if a solution curve hits a nullcline somewhere other than an equilibria it immediately crosses to the other side of the nullcline, tangency is not possible since the nullclines are nowhere vertical (for the one representing vertical solution curves) and nowhere horizontal (for the other).

The inside region
Consider the region formed by the segment of the parabola P 1 from S to the first intersection with the line w 2 = 1 (this is at T if θq1 r1 > 1); if the intersection is not at T then use the line w 2 = 1 to connect the endpoint of the segment to T ; the segment of P 2 from S to the first intersection with the line w 1 = 1 (this is at T if θq2 r2 > 1) and if the intersection is not at T then use the line w 1 = 1 to connect the endpoint to T . Thus, the inside region is defined to be the region inside both parabolae and the unit square, in the notation of the previous section it is a subset of ω 2 . See Figures 4 and 5.
Denote the edges of the region as follows: • The open segment of P 1 (not including S or the other endpoint) by A.
• The open segment of P 2 (not including S or the other endpoint) by B.
• If A does not connect S to T then denote the open segment of w 2 = 1 (not including endpoints) by C. Notice that a monotone connection from S to T must approach T from inside this region if both C and D exist. If either C or D does not exist then it is possible to reach T monotonically from outside the inside region -then we must use the outside region as discussed below.
If C exists, then denote its left endpoint by x and if D exists then denote its lower endpoint by y.
By (20) and the following discussion, a necessary condition for a monotone connection that eventually approaches T from the inside region is that since w 1 > 0 and w 2 > 0 (throughout the interior of the region) if and only if (C1) holds.

Lemma 4.1 Assume (C1) holds. Then a solution curve that hits the boundary of the inside region at any point other than S or T immediately crosses the boundary.
Proof. Consider the various segments of the boundary. On A and at the point x, if it exists, w 1 = 0 and w 2 > 0. Thus (using continuity and the fact that the slope of this boundary is bounded and positive, so the flow is not tangent to the boundary) if w(t) ∈ A or w(t) = x then there exists an > 0 such that for s ∈ (t − , t), w(s) is inside the region and for s ∈ (t, t + ), w(s) is outside the region.
Similarly for B and y -w 2 = 0 and w 1 > 0 so the same argument works.
For C, if it exists, note that w 1 > 0 and w 2 > 0 on C and again the solution curve crosses from inside to outside automatically. Similarly for D, and the lemma is proved.

The outside region
We have observed that there are monotone curves from S to T that do not intersect the inside region, except at T itself, when the inside region does not have four edges. Without loss of generality assume that the edge C of the inside region does not exist, i.e. θq1 r1 > 1, so that there is a possibility of a monotone connection through ω 3 . When neither C nor D exists there are two outside regions, by equation (20) the flow is in opposite directions in the two regions. Hence there can be a monotone connection through at most one of them.
Consider the region formed by the segment of the parabola P 1 from S to T ; the w 2 axis and the line w 2 = m(w 1 − 1) + 1 for some m ∈ (0, 1 − r1 θq1 ). This restriction on m ensures that the line has positive slope and lies above the parabola P 1 (the upper bound just given is the slope of the parabola at T and the assumption that θq1 r1 > 1 ensures that the slope is strictly positive). We will choose a particular value of m later. In the notation of section 4.2 this region is therefore a subset of ω 3 . See Figure 6 for a diagram of a typical example of this region.
Denote the edges of the region as follows: • The open segment of P 1 (not including S or T ) by A. • The open segment of w 2 = m(w 1 − 1) + 1 strictly between the w 2 axis and T by C.
For a monotone connection from S to T that does not go through the inside region it is necessary that w 1 > 0 and w 2 > 0 throughout the interior of the outside region (note that though there are curves from  Figure 6: An example of an outside region S to T that exit the outside region through C and still converge monotonically to T this condition will still be necessary).
Using equation (20) to consider the direction of the flow within the region notice that is is therefore necessary that b 2 > b 1 and that −b 2 < c < −b 1 , i.e. (b 2 + c) > 0 > (b 1 + c). So, for there to be a possibility of a monotone connection that is not through the inside region, the following conditions are necessary: and For a monotone connection through the other outside region it is necessary that (b 1 +c) > 0 > (b 2 +c) and θq2 r2 > 1 -the subsequent analysis is entirely equivalent since we can simply interchange the subscripts 1 and 2. The two cases are clearly mutually exclusive.
The condition (C2) ensures that, for each k > 0, a curve that satisfies w 2 /w 1 = k is an ellipse passing through S and T . This observation is used below to rule out internal tangencies to the region.
A solution curve cannot exit from the outside region through A or B in forwards time. On A observe that w 1 = 0 and w 2 > 0 and on B that w 1 > 0 and w 2 > 0.
We need to check that the flow does not have an internal tangency to C in order to construct a shooting argument for the outside region. The shooting argument used runs backwards in time, but note that a tangency backwards in time is also a tangency forwards in time, and vice versa.

Lemma 4.2 Assume condition (C2) holds. Then a solution curve that hits the boundary of the outside region (from the inside of the region -an external tangency is immaterial) at any point other than S or T immediately crosses the boundary.
Proof. If a solution curve hits A or B the observations above show that the curve (in backwards time) will cross from the inside to the outside of the region. So assume, for a contradiction, that at a point x 1 on C the solution curve is internally tangent to C. The solution curve in forwards time will re-enter the interior of the region -where w 1 > 0 and w 2 > 0 -and must eventually hit the line C again or else approach T . This is since it cannot leave the region through A because w 1 = 0, w 2 > 0 there. Call the first point at which this occurs x 2 .
Thus, the solution curve is a smooth curve from x 1 to x 2 , both of which are on the line w 2 = m(w 1 −1)+1, therefore there is a point x 3 on this curve at which the curve has slope m. However, the curve on which the slope of the flow is m is an ellipse; an ellipse that passes through S, T and x 1 -so x 3 cannot be on this ellipse. This contradiction proves the result.

The shooting argument
To establish the main result (Theorem 1.1) we shall use the preceding observations in a shooting argument through the outside region (backwards in time), as well as a more standard forward time argument from S to T for the inside region. Note that a monotone travelling wave from S to T is possible only if T has stable monotone eigenvalues. A corresponding unstable direction at S is also necessary, but this always exists when a stable monotone eigenvalue at T exists. (See the discussion before Lemma 3.2 in section 3.) It is known so far that if c > max(−b 1 , −b 2 ) then a monotone connection can only exist through the inside region -we show below that in this case it does indeed exist and is unique. For c < min(−b 1 , −b 2 ) there cannot be a monotone connection from S to T and for c in between then there may or may not be a monotone connection (which necessarily lies entirely in an outside region) -the determining factor is whether c is greater than c(θ). When such a connection exists it is unique.

Through the inside region when (C1) holds
Consider the intersection of the inside region with the circle of radius , centred at the origin, S. This is a segment of a circle with 2 endpoints on the boundary of the inside region. w 1 > 0 and w 2 > 0 are necessary for a monotone connection through the region (i.e. (b 1 + c) > 0 and (b 2 + c) > 0 are necessary in this region) -when this is true a connection looks plausible. Considering the family of solution curves which pass through points of this segment at time zero completes our shooting argument, since, running backwards in time all these curves must go to S (since they cannot exit the region and cannot tend to any point other than an equilibrium), and running forwards in time, the curves from the 2 endpoints leave the region immediately and curves from points in between will exit transversally from the region (by Lemma 4.1), thus exiting in between the exit/end-points. Hence classical continuous dependence theory for initial value problems allows us the conclude that a monotone connection from S to T exists.

Through the outside region when (C2) holds
Assume condition (C2), otherwise there is no monotone connection through the region. This time we shoot backwards in time.
Firstly we deal with c > c(θ). Choose m so that the edge C has slope equal to half that of the dominant eigenvector of K c,θ (T ) (the dominant eigenvector is that corresponding to the negative real eigenvalue of smallest modulus -when (C2) and c > c(θ) there are two simple stable monotone eigenvalues so m is well-defined). Thus the upper edge of the region bisects the angle between the line w 2 = 1 and the dominant eigenvector at T . This enables us, backwards in time, to obtain paths leaving the outside region on both sides of T , and hence shoot to S. More precisely, for c > c(θ) we already know that K c,θ (T ) has two stable monotone eigenvalues. Consideration of equation (17) shows that both have slope less than that of the P 1 at T . The dominant eigenvector is 1 v + . Thus we choose m to correspond to a direction 1 v + /2 . Hence, by shooting backwards from points suitably close to T and using the fact that the flow is determined by the dominant eigenvector we observe the flow sweeps through the boundaries of the region (by Lemma 4.2) as we follow a segment of a circle around T that intersects the region. Classical continuous dependence theory for initial value problems tells us that a connection from S to T exists.
For c < c(θ) there is no stable monotone eigenvalue at T and there cannot be a monotone connection. Thus we have completed showing that a monotone eigenvalue at T is sufficient (except for the cases c = c(θ) and c = max(−b 1 , −b 2 ) > c(θ) which are discussed below in section 4.6) for a monotone connection -we already knew it was necessary.

Special cases
We have been wary of the cases where c = c(θ) and c = −b i . Assume, without loss of generality, that b 1 > b 2 .
If θq 2 ≤ r 2 , then c(θ) = −b 2 . For c = c(θ) the only candidate for a connection from S to T is the segment of the w 2 -nullcline, P 2 , connecting S to T . If θq 2 < r 2 then this is not monotone from S to T and so there is no monotone travelling wave at this speed. When θq 2 = r 2 this connection is monotone. The travelling wave equations consists of one algebraic -defining the nullcline -and one differential equation. We can substitute the algebraic into the differential equation to obtain a one-dimensional problem where we are looking for a connection from 0 to 1. Since the derivative is positive between these points this indeed exists and the segment of nullcline is a monotone travelling wave. When For the case c = c(θ) we can use the argument we used in the outside region -by part (iii) of Lemma 3.2 there is a double stable monotone eigenvalue so the corresponding eigenvector determines the nature of the flow near T . We simply repeat the bisection procedure and show that the flow leaves either side of the region -so there is a monotone connection for c = c(θ) in this case. We also must check the case c = −b 2 , but here the solution is simply the segment of nullcline exactly as for θq 2 = r 2 .
Thus we can summarize as follows: • For θq 2 < r 2 , there is a monotone travelling wave if and only if c > c(θ).
• For θq 2 ≥ r 2 . There is a monotone travelling wave if and only if c ≥ c(θ).

Uniqueness modulo translation
To establish uniqueness of monotone travelling waves from S to E i , i = 1, 2, 3, it suffices, because of the change of variables at (1.26), to prove that the wave from S to T is unique. However, due to the dimensionality of the unstable and stable manifolds at S and T respectively, we have very little left to show. As we summarised in section 3; just above Lemma 3.2; we have that: • For c < min(−b 1 , −b 2 ); K c,θ (S) has a two-dimensional stable manifold. (Of course, here c < c(θ) anyway.) ; K c,θ (S) has a one-dimensional unstable manifold.
Thus, since we can only possibly have non-uniqueness when: • K c,θ (S) has a two-dimensional unstable manifold and K c,θ (T ) has a two-dimensional stable manifold, we have only one case left to check, i.e. c > max(−b 1 , −b 2 ) and θ < 1 ρ1+ρ2 . In this case the two eigenvalues of T are negative and distinct, denote them by α < β < 0, and the corresponding eigenvectors by v α and v β , respectively. α is a stable monotone eigenvalue (i.e. v α has both components of the same sign) and β is not. Note that the dominant eigenvalue β is the nonmonotone eigenvalue. This suggests that most of the flow entering T will do so non-monotonically and so there may indeed be a unique monotone connection. We now use this idea to complete the proof of uniqueness.
Firstly, by standard theory (Coddington and Levinson [3,Chapter 13,Theorem 4.4]), all solutions that converge to T do so exponentially. For notational convenience define u by w i (x) = u i (x)+1 for i = {1, 2}, so that u converges to the origin when w converges to T . Then there exists a one-dimensional manifold of solutions u such that and all other solutions that converge to the origin do so with rate β. Hence uniqueness is proven and we have completed the proof of Theorem 1.1.

Martingales for a branching two-type process
Consider the two-type branching system of particles which was defined in section 1.4 and recall our notation that I = {1, 2} and that we often switch between two equivalent notations in this section to reduce the use of subscripts. Thus we will sometimes write, for example, b(y) for b y , w(x, y) for w y (x) and u(t, x, y) for u y (t, x) (for y = 1, 2).
The formal generator G of the two-type branching process is where, for n ≥ 1, x ∈ R n , and y ∈ I n , we have  (7)) and ∂ ∂t + G F (t; n; x; y) = 0 (n ≥ 1, x ∈ R n , y ∈ I n ), then F t; N (t); X(t); Y (t) is a local martingale. In particular, if u is a C 1,1 solution of the coupled system (1), define, for t > 0, This satisfies the conditions of equation (21) (the partial differential equation enables us to do the time differentiation of M ) and hence defines a 'multiplicative' local martingale M on the time-parameter set [0, t]. If in addition 0 ≤ u(·, x, y) ≤ 1, then 0 ≤ M ≤ 1, so that M is a true martingale and that is, Thus we have proved the following: then u has a McKean representation The above theorem clearly tells us that there is at most 1 bounded (between 0 and 1) solution to the system (1) for smooth, bounded (between 0 and 1) initial data. Now we look at 'additive' martingales.
Now we have: is a local martingale if, and only if, w solves the travelling-wave system (2). If g is a C 1 function on R × I, then Then, by Theorem 5.2, Z λ is an 'additive' local martingale; and since it is non-negative, it is also a supermartingale. We wish to show that Z λ is a true martingale, this can be done just as for ζ λ by noting that EN (t) ≤ e r0t and each individual term of the sum is again bounded up to any fixed t.

A large-deviations approach to c(θ)
As we have already observed in Theorem 1.3, the asymptotic behaviour of the position of the left-most particle of the system should give us the wave speed c(θ). This is proved later in this section; for now we want to give probabilistic heuristics to obtain this value.
If we consider a particle's type intuitively we note that it flicks between the two types with an equilibrium distribution independent of θ. On average the particle will spend a proportion q2 q1+q2 of its time with type 1 and a proportion q1 q1+q2 of its time with type 2. A larger θ will tend to force each particle's type distribution closer to the equilibrium distribution.
So, we are interested in considering how far from the equilibrium distribution (for time spent in the two states) a particle can deviate. The largest deviations will give the positions of the most extreme particles, including the left-most particle. Clearly, the higher the breeding rate (the parameters r 1 , r 2 ) for each state, the easier it will be to find a particle in that state. As θ is increased the deviations from the equilibrium distribution will be smaller.
The theory of large deviations is the appropriate technique for tackling such a problem. We use the Dirichlet form (u, v) := −u T ΠQv, where Π is the diagonal matrix whose elements are the components of the equilibrium distribution for the Markov chain. This is a positive-definite, symmetric bilinear form. Then the large-deviation rate functional J (µ) is (f 1/2 , f 1/2 ) where f is the ratio of µ to the invariant distribution and µ is any distribution for time spent in the two states (see Deuschel and Stroock [7,page 129]). So, if µ corresponds to spending a proportion of time m 1 in state 1, then J (µ) = ( √ q 1 m 1 − √ q 2 m 2 ) 2 , and hence enables us to state that the probability that a particle, at time t, has spent a proportion m 1 of its time in state 1 and a proportion m 2 of its time in state 2 (for m 1 , m 2 ≥ 0, Thus, for all proportions other than the equilibrium proportions this decays exponentially and an increase in θ increases the rate of decay. Now, this decay is balanced by the breeding -the number of particles grows exponentially at a rate m 1 r 1 + m 2 r 2 . So, when there is exponential growth in the number of particles who have spent a proportion m 1 of their time in state 1 and a proportion m 2 of their time in state 2; when this is strictly less than zero there is decay in the number of such particles. This leads us to definec(θ) as: where m 1 and m 2 represent the proportions of times in the two states. This is well-defined since for m 1 = q2 q1+q2 (the equilibrium distribution) the inequality is certainly satisfied. To demonstrate the correspondence of c(θ) andc(θ) we examine some specific cases.
Firstly if b 1 = b 2 = b then the above formula givesc(θ) = −b, which tallies with the result for c(θ) we obtained in section 3. As before we can now assume, without loss of generality, that b 1 > b 2 . For the equilibrium proportions the inequality (23) is satisfied, so we immediately have that So we ask, when isc(θ) = −b 2 ? Putting m 2 = 1 the inequality (23) becomes r 2 − θq 2 > 0, so r 2 > θq 2 implies thatc(θ) = −b 2 . For r 2 = θq 2 we use the concavity of m 1 r 1 + m 2 r 2 − θ( √ q 1 m 1 − √ q 2 m 2 ) 2 to observe that the inequality (23) is satisfied for m 2 = 1 − for all in the interval 0 < < q2 q1+q2 , thus the infimum is still m 2 = 1 -just it is not attained (this corresponds to the fact that in this case there is a travelling wave corresponding toc(θ) = c(θ) , for r 2 > θq 2 there is not).
For r 2 < θq 2 ,c(θ) < −b 2 and so all that remains to prove is thatc(θ) = c(θ) in this case too. We can rearrange our formula forc(θ) by writing c = −b 1 m 1 + b 2 (m 1 − 1) and hence m 1 = b2+c b2−b1 . Substituting this into the defining inequality we note that implies that the term under the square root in equation (16) is zero. Now we note that this corresponds to the zero with the larger value of c sincec(θ) is at least c * .
This formula for c(θ) makes it easier to observe that we can have m i = 1 only when r i > θq i , which makes intuitive sense -the breeding in state i needs to be larger than the mutation out of the state for there to be a persistent family of that type of particle over time. It also makes it clear that c(θ) decreases from max(−b 1 , −b 2 ) to c * as θ → ∞.

Convergence properties of Z λ martingales
The full assumptions which we have made about Z λ (defined in Theorem 1.4) -that c > c(θ), that λ is the probabilistic eigenvalue of K c,θ , etc. -will now be needed in proving that Z λ converges in L p for some p > 1 (and hence in L 1 ). According to Doob's L p inequality (see Rogers and Williams [14], Theorem II.70.2), we need only show that Z λ is bounded in L p for some p > 1. The key result is from Neveu [13] and the method of using it follows Champneys et al. [2].
Proof of L 1 convergence of Z λ . Fix t > 0. Because of the branching character of the (N, X, Y ) process, we have for each s > 0, where, conditionally on F s , the W k (t, s) are independent, each with the P 0,y(k) law of Z λ (t) where y(k) = Y k (s) ∈ I. Since t is fixed, and I is finite, Neveu's lemma, applied conditionally on F s , gives so that, on taking expectations, where µ := λp < λ < 0. We choose p in (1,2] sufficiently close to 1 that (see Lemma 3.2) c > c 1 , where −µc 1 = Λ P F (µB + θQ + R) and v µ is the corresponding eigenvector with v µ (1) = 1. But then and Z λ is bounded in L p . Now we prove, with the same notation and assumptions, that w(y) := P x,y Z λ (∞) = 0 = 0 for all (x, y).
(The fact that w(y) does not depend on x is obvious.) Proof of (24). Let J be the first jump time of Y 1 and let T be the first branch time of (N, X, Y ). On decomposing w(y) according as T < J or T > J, we obtain w(y) = r(y)w(y) 2 + θ z =y Q(y, z)w(z) r(y) + θq(y) , so that Rw = R(w 2 ) + θQw. By Lemma 3.4, w ≡ 0 on I or w ≡ 1 on I. When Z λ converges in L 1 , then, obviously, w ≡ 0 on I.
The following theorem is now proven.

Theorem 5.4 Let c > c(θ).
Let λ be the probabilistic eigenvalue of K c,θ (T ). Thus −λc is the Perron-Frobenius eigenvalue of λB + θQ + R and, as usual, we denote by v λ the corresponding eigenvector with v λ (1) = 1. Then is a true (not just a local) martingale, and Z λ (t) converges to a limit Z λ (∞) almost surely and in L 1 . Moreover, P x,y (Z λ (∞) > 0) = 1 for all x and y.
For the proof of the uniqueness modulo translation of the monotone travelling wave from S to T , we need the following result. When there is only one stable monotone eigenvalue of K c,θ (T ) we have no other candidate for a wave (each stable monotone eigenvalue gives us a martingale which could correspond to a suitable wave), when there are two we use this result to rule out the other possibility.

Lemma 5.5
Suppose that c > c(θ) and that there are two stable monotone eigenvalues of K c,θ (T ). Denote by β the eigenvalue further from 0, and by v β the associated Perron-Frobenius eigenvector of βB + θQ + R with v β (1) = 1. Then, almost surely, Proof. We prove this result by modifying an argument in Neveu [13]. Let 0 < p < 1. Then for u, v > 0, Again, let J be the first jump time of Y 1 and let T be the first branch time of (N, X, Y ). The decomposition leads to the formula where α := pβ. On evaluating these expectations, we obtain g(y) ≤ θ z =y Q(y, z)g(z) + 2r(y)g(y) − b(y) + c α + r(y) + θq(y) , which rearranges to give 0 ≤ (B + cI)α + θQ + R g.
We know that g ≥ 0 on I. Lemma 3.3 shows that if we choose p sufficiently close to 1, then g = 0. The lemma is proved.

Proof of Theorem 1.4
We now verify steps in the proof of Theorem 1.4, which lead to the important equality (8) of Theorem 1.3: Part (i). Let c > c(θ), and let λ be the probabilistic eigenvalue of K c,θ (T ). Let Z λ be the associated martingale. We see by considering the position of the left-most particle that for finite constants K 1 (λ), K 2 (λ) and K 3 (λ), independent of t. Hence, for > 0, by Doob's submartingale inequality. By the Borel-Cantelli lemma, we have t −1 M λ (t) → 0, a.s., whence, from (25), lim sup Part (ii) of Lemma 3.2 clinches Part (iii) of Theorem 1.4, and the proof of (8) is complete.

Heaviside initial conditions
We have the McKean representation for the unique (by the work of section 2) solution of our coupled equation (1) when 0 ≤ u ≤ 1 and the initial data u(0, ·, ·) are sufficiently smooth. We would like to obtain a McKean representation in the case of the Heaviside initial data. We verify this directly below then, because t −1 L(t) → −c(θ) (a.s.), the rest of Theorem 1.3 is obvious. In a separate paper (currently in preparation), Lyne and Williams will • give a new direct analytic method of proving existence and uniqueness of weak solutions with bounded, measurable initial data, • show that these solutions lead to martingales and hence to a McKean representation, • give numerical and simulation studies of the solutions with Heaviside initial data, and • study the long-term behaviour of these Heaviside solutions.
We can observe several features of the solution from this representation. Since the left-most particle must lie in the interval [x + t min(b 1 , b 2 ), x + t max(b 1 , b 2 )] at time t then the solution is identically zero for x ≤ −t max(b 1 , b 2 ) and is identically one for x > −t min(b 1 , b 2 ). It is clear that the solution retains the left-continuity of the initial data. Also, u(t, x, y) must be strictly increasing (in x) inside this region.
The case b 1 = b 2 = b, say, is very simple. All particles will travel at speed b, so that, at time t, all particles will be at x + bt. Hence, L(t) = x + bt, so P x,y [L(t) > 0] = I {x>−bt} . This is simply the Heaviside function travelling at speed −b. This corresponds to the analysis in section 2.6.
To deal with b 1 = b 2 assume from now on, without loss of generality, that b 1 > b 2 .
Note that the only way that the left-most particle can be at position x + b 1 t at time t > 0 is if the first particle was of type 1 and there have been no mutations to type 2 by the particle or any of its descendants so far. The probability of this event occurring when the initial particle is of type 1 is which explains why u(t, −t+) equals this in the calculation in section 2.4. The only way that the left-most particle can be at position x + b 2 t at time t > 0 is if the first particle was of type 2 and the embedded birth and death process on this line has not become extinct. This embedded process is that constructed by considering a particle to die when it mutates to type 2. Thus, on x = −b 2 t, for t > 0, u 1 is identically zero. The probability of this event occurring when the initial particle is of type 2 is exactly the expression calculated for v(t, +t) in section 2.4, and the limit as t → ∞ (which is min( θq2 r2 , 1)) is the probability of eventual extinction. More generally it is natural to expect that the McKean representation will tie up with weak solutions of our coupled equation (1). This is because they are piecewise classical and any discontinuities (corresponding to atoms of probability in the distribution of the particles' positions) will naturally propagate only along characteristics (since the atoms of probability travel at b 1 and b 2 ). Thus a solution satisfying the McKean representation satisfies the Rankine-Hugoniot conditions, hence is a weak solution (see section 2.3). This function w satisfies the travelling-wave equation (2) and is, modulo translations, the unique monotone wave of speed c from S to T .

Polishing off the probability
Proof. We are guided by McKean [11]. So, we are supposing that u solves (1), that 0 ≤ u ≤ 1 and that 1 − u(0, r, y) ∼ v λ (y)e λr (r → ∞). Since Z λ (∞) exists almost surely (and clearly all the terms in the sequence of inequalities are bounded below by 0 and above by 1) this becomes After taking expectations and using Fatou's Lemma, this implies that On letting ↓ 0, we now obtain the desired result u(t, x + ct, y) → w(x, y) = E x,y exp{−Z λ (∞)}.
Existence of a monotone travelling wave from S to T when c > c(θ). It is now intuitively obvious, and not that difficult to prove directly from the branching property, that the function w(·, ·) in (28) is a monotone travelling wave from S to T . Firstly it does in fact satisfy the travelling wave equation (2)  w(x, y) = r(y) r(y) + θq(y) w 2 (x, y) + θq(y) r(y) + θq (y) w(x,ŷ) + b(y) + c r(y) + θq (y) w (x, y) since the last term in the integration by parts is simply a multiple of the x-derivative of equation (29).
Hence our claim is proved.
Uniqueness modulo translation of the monotone travelling wave from S to T . Let c > c(θ), and letw be a monotone travelling wave from S to T . We know from differential-equation theory that either a suitable translate ofw satisfies (30), or a suitable translate ofw satisfies where β is the monotone eigenvalue of K c,θ further from 0.
The proof of uniqueness is now complete, so we have proven Theorem 5.6.
It is interesting to compare the above probabilistic proofs of existence and uniqueness modulo translation of travelling waves with the analytic proofs given in section 4. We could, for example, use ODE results to obtain results on L 1 convergence of our martingales.
The methods in section 4 dealt with the c = c(θ) case; the probability theory for this case should be amenable to the techniques of Neveu [13] -stopping lines will correspond to simple conditions on the occupation times of the two states of the Markov chain.