On Dirichlet eigenvectors for neutral two-dimensional Markov chains

We consider a general class of discrete, two-dimensional Markov chains modeling the dynamics of a population with two types, without mutation or immigration, and neutral in the sense that type has no influence on each individual's birth or death parameters. We prove that all the eigenvectors of the corresponding transition matrix or infinitesimal generator \Pi\ can be expressed as the product of"universal"polynomials of two variables, depending on each type's size but not on the specific transitions of the dynamics, and functions depending only on the total population size. These eigenvectors appear to be Dirichlet eigenvectors for \Pi\ on the complement of triangular subdomains, and as a consequence the corresponding eigenvalues are ordered in a specific way. As an application, we study the quasistationary behavior of finite, nearly neutral, two-dimensional Markov chains, absorbed in the sense that 0 is an absorbing state for each component of the process.


Introduction
This paper studies spectral properties of two-dimensional discrete Markov processes in continuous and discrete time, having the neutrality property, in the sense of population genetics (see e.g. [17]). Considering two populations in interaction, corresponding to two different types of individuals (typically a mutant type and a resident type), one says that the types are neutral (or the mutant type is neutral w.r.t. the resident type, or more simply the mutation is neutral) if individuals of both types are indistinguishable in terms of the total population dynamics. In other words, the mutant population has no selective advantage (or disadvantage) with respect to the rest of the population.
More formally, we say that a two-dimensional Markov process (X n , Y n ) n∈T (where T = Z + ≔ {0, 1, . . .} or R + ) with values in R 2 + or Z 2 + is neutral if Z n = X n + Y n is a Markov process. In particular, the law of the process Z depends on Z 0 , but not on X 0 or Y 0 .
If the process Z is a birth and death continuous-time chain, the class of neutral processes we consider is the following: the birth and death rates of the Markov process (Z t ) t∈R + when Z is in state k ≥ 0 are of the form kλ k and kµ k , respectively. Note that 0 is an absorbing state for Z. With this notation, the parameters λ k and µ k can be interpreted as birth and death rates per individual. Then the process (X t , Y t ) t∈R + is the birth and death process where both types of individuals have birth and death rates per individual λ k and µ k , when the total population has k individuals. This leads to the following transition rates for the Markov process (X t , Y t ) t∈R + : for all (i, j) ∈ Z 2 + , from (i, j) to (i + 1, j) with rate i λ i+j from (i, j) to (i, j + 1) with rate j λ i+j from (i, j) to (i − 1, j) with rate i µ i+j from (i, j) to (i, j − 1) with rate j µ i+j .
Note that the sets {0} × Z + , Z + × {0} and {(0, 0)} are absorbing for this process. In other words, we only consider neutral two-dimensional processes without mutation and immigration.
In the case of discrete time, we consider two-dimensional birth and death processes constructed in a similar way: assume that the birth and death probabilities of the process (Z n ) n∈Z + when in state k are p k and q k , respectively, with p k + q k ≤ 1. Then, when a birth or a death occurs in the population, the individual to which this event applies is chosen uniformly at random in the population. This leads to the transition probabilities from (i, j) to (i + 1, j) with probability i i+j p i+j from (i, j) to (i, j + 1) with probability j i+j p i+j from (i, j) to (i − 1, j) with probability i i+j q i+j from (i, j) to (i, j − 1) with probability j i+j q i+j from (i, j) to (i, j) with probability r k , where r k ≔ 1 − p k − q k . Note that this construction requires assuming that r 0 = 1 (i.e. that 0 is absorbing for Z).
In [14], Karlin and McGregor studied two families of neutral multitype population processes (branching processes and Moran model), but only in the case of nonzero mutation or immigration, for which the set of states where one population (or more) is extinct is not absorbing. They could express the eigenvectors of the corresponding infinitesimal generators in terms of Hahn polynomials. We focus here on neutral processes without mutation and immigration, which are singular for the approach of [14], and we apply our study to a much bigger class of neutral population processes, containing the birth and death processes described above, but also non-birth and death models.
Our main result is the characterization of all eigenvalues and right eigenvectors of the transition matrix of neutral processes without mutation and immigration. To this aim, we first consider the (easier) continuous state space case in Section 2 to introduce some tools used in the sequel. Next, we construct a particular family of polynomials of two variables in Section 3, using linear algebra arguments. In Section 4, we prove that the eigenvectors of the transition matrix of neutral two-dimensional Markov processes can be decomposed as the product of "universal" polynomials (in the sense that they do not depend on the specific transition rates of the Markov chain) with functions depending only on the total population size. We then relate these eigenvectors with Dirichlet eigenvalue problems in subdomains of Z 2 + of the form {(i, j) ∈ N 2 : i + j ≥ k} for k ≥ 2 (Section 5), where N = {1, 2, . . .}.
The last section (Section 6) is devoted to the application of the previous results to the study of quasi-stationary distributions. A probability distribution ν on Z 2 + \ {0} is called quasi-stationary if it is invariant conditionally on the non-extinction of the whole population, i.e. if P ν ((X 1 , Y 1 ) = (i, j) | Z 1 = 0) = ν i,j , ∀(i, j) ∈ Z 2 where P ν denotes the law of the process (X, Y ) with initial distribution ν. This question is related to the notion of quasi-limiting distribution (also called "Yaglom limit", in reference to Yaglom's theorem on the same convergence for Galton-Watson processes), defined as These notions are relevant in cases where extinction occurs almost surely in finite time, to describe the "stationary behaviour" of the process before extinction when the extinction time is large. This is typically the case in many population dynamics models, where ecological interactions in the population produce high mortality only when the population size is large (one speaks of density-dependent models, see e.g. [18], or [5] for discrete stochastic models). These questions have been extensively studied in the case where the transition matrix restricted to the non-extinct states is irreducible (which is not true in our two-dimensional case). The first paper of Darroch and Seneta [6] studies the discrete-time, finite case. Several extensions of these results to continuous-time and/or infinite denumerable state spaces have then been considered in [21,7,10]. The case of population dynamics in dimension 1 have been studied by many authors(e.g. [4,22,16,9,15,13,19,12]).More recently, the quasi-stationary behaviour of one-dimensional diffusion models has been studied in [2]. As far as we know, the two-dimensional case has only been studied in the continuous state space (diffusion) case [3]. An extensive bibliography on quasi-stationary distributions can be found in [20] In Subsection 6.1, we first give the quasi-limiting distribution for general finite two-dimensional Markov chains in terms of the maximal Dirichlet eigenvalues of the transition matrix in several subdomains. Finally, in Subsection 6.2, we apply our previous results to prove that coexistence in the quasi-limiting distribution is impossible for two-dimensional finite Markov chains which are close to neutrality.
The paper ends with a glossary of all the notation used in Sections 4 to 6, which may appear at different places in the paper.

Preliminary: continuous case
In this section, we consider the continuous state space case, where computations are easier, in order to introduce some of the tools needed in the discrete case.
Fix p and q two measurable functions from R + to R + , and consider the system of stochastic differential equations where (B 1 , B 2 ) is a standard two-dimensional Brownian motion. Such SDEs are sometimes called branching diffusions, and are biologically relevant extensions of the classical Feller diffusion [2]. If p and q satisfy appropriate growth and regularity assumptions, the solution to this system of SDEs is defined for all positive times and can be obtained as scaling limits of twodimensional birth and death processes (we refer to [2] for the one-dimensional case; the extension to higher dimensions is easy). This process is neutral in the sense defined in the introduction since where is a standard Brownian motion. Note also that R + × {0}, {0} × R + and {(0, 0)} are absorbing states as soon as there is uniqueness in law for the system (1). For any ϕ ∈ C 2 (R 2 + ), the infinitesimal generator A of the process (X t , Y t ) t≥0 is given by We first observe in the following proposition that A admits a symmetric measure, but only on a subset of C 2 (R 2 + ).
Proposition 2.1 Assume that 1/p and q/p belong to L 1 loc ((0, +∞)). Let us define the measure µ on R 2 + as Then, the restriction A of the operator A to C 2 c ((0, +∞) 2 ) is symmetric for the canonical inner product ·, · µ in L 2 (R 2 + , µ), and, hence, so is its closure in L 2 (R 2 + , µ). Note that, because of the singularities in µ when x or y vanish, if p ≥ c > 0 in the neighborhood of 0, any continuous function in L 2 (R 2 + , µ) must vanish at the boundary of R 2 + . Therefore, Proof For all f, g ∈ C 2 c ((0, +∞) 2 ), we have (formally) Conversely, these equalities can be directly checked from the formula (4), which implies that ∂(xµp)/∂x and ∂(yµp)/∂y exist in a weak sense.
Before studying the eigenvectors of A, we need the following result.
has no (weak) non-zero solution h ∈ C 1 ([−1, 1]) except when λ = d(d − 1) for some d ∈ N. For d = 1, the vector space of solutions has dimension 2 and is spanned by the two polynomials h(x) = 1 and h(x) = x. For all d ≥ 2, has a one-dimensional vector space of solutions in C 1 ([−1, 1]), spanned by a polynomial H d of degree d, which can be chosen such that the family (H d ) d≥2 is an orthonormal basis of L 2 ([−1, 1], dx 1−x 2 ). In addition, H d has the same parity as d (all powers of x that appear are even if d is even and odd if d is odd).
Proof The case λ = 0 is trivial, so let us assume that λ = 0. Note that, assuming λ = 0 implies that any solution h ∈ C 1 ([−1, 1]) satisfying (5) must satisfy h(−1) = h(1) = 0. Indeed, the function h(x)/(1 − x 2 ) has a finite integral over [−1, 1] in that case only. The equation has the form of the classical Sturm-Liouville problem, but does not satisfy the usual integrability (5) can be proved using the classical Wronskian determinant: given two solutions h 1 and h 2 of (5) in ) are linearly dependent, and so h 1 and h 2 as well.
Second, assuming that (5) has a polynomial solution of degree d, the comparison of the higher degree terms yields λ = d(d − 1). So assume that where a d,d = 0 and is a nonzero solution of (6). Since 1], and for any d ≥ 2, where the integration by parts is justified by the Dirichlet conditions for h and H d . Therefore, either λ = d(d − 1) for some d ≥ 2, and then f and H d are linearly dependent by the previous uniqueness result, or is a family of polynomials of all nonnegative degrees. Therefore, 1 −1 h(x)P (x)dx = 0 for all polynomial P , and hence h = 0.

Remark 2.3
The polynomials H d may be obtained as a limit case of Gegenbauer polynomials, defined for a parameter λ > −1/2 as (up to a multiplicative constant) for all d ≥ 0, where B is the Beta function, a F b are the classical hypergeometric series and (x) n := x(x + 1) . . . (x + n − 1) is the shifted factorial. The normalizing constant has been chosen such that these polynomials form an orthonormal family in L 2 ((−1, 1), (1 − x 2 ) λ−1/2 dx). Then, for all d ≥ 2, H d may be obtained as One easily checks from this formula that H d solves (6), and the fact that We now introduce the change of variables which defines a C ∞ -diffeomorphism from R 2 Now, assume that A has an eigenvector of the form ϕ(w)ψ(z). If λ is the corresponding eigenvalue, (11) yields Hence, we have the following result.
Proposition 2.4 All functions on R 2 + of the form where for some λ ∈ R, are eigenfunctions of A for the eigenvalue λ. Now, we proceed to prove that there are no other eigenvectors of A (in the sense of Theorem 2.5 below). Such a result seems natural as (13) can be written as the Sturm-Liouville problem Here again, the usual integrability conditions are not satisfied. More precisely, if p(z) ≥ c > 0 for z in the neighborhood of 0 and q is bounded, using the terminology of [23], for all d ≥ 0, the problem (12) is a singular self-adjoint boundary value problem on (0, ∞), where the endpoint 0 is LP (limit point) singular (see [23,Thm. 7.4.1]). In this case very little is known about the existence of an orthonormal basis of eigenvectors in L 2 ((0, ∞), ν), where (the spectrum might even be continuous, see [23,Thm. 10.12.1(8)]). For this reason, we state our next result on the operator A with a restricted domain corresponding to the case where the diffusion is reflected in the set t . Then Z t = X t + Y t is the solution of (2) reflected at a and b with local time The corresponding infinitesimal generator is defined by (3) with domain the set of ϕ(x, y) ∈ C 1 (D) ∩ C 2 (int(D)), where int(D) denotes the interior of D, such that Theorem 2.5 For 0 < a < b < ∞, assume that p ≥ c > 0 on [a, b] and q/p belong to L 1 ([a, b]).
(a) There exists a denumerable orthonormal basis of L 2 (D, µ) of eigenvectors of A of the form (12), where d ≥ 2 and ψ solves (13) on (a, b) and satisfies ψ ′ (a) = ψ ′ (b) = 0. Moreover, any eigenvector of A in L 2 (D, µ) is a linear combination of eigenvectors of the form (12), all corresponding to the same eigenvalue.
(b) There exists a family of right eigenvectors of A of the form which is a basis of the vector space where Leb denotes Lebesgue's measure. More precisely, for all f ∈ V , the functions f 1 , f 2 , f 3 in (16) are unique and there exists unique sequences {α k } k≥1 , {β k } k≥1 , and {γ dk } d≥2, k≥1 such that where the series k α k ψ k and k β k ψ k both converge for · ν and d, k (x + y) converges for · µ . Point (b) says that the eigenvectors of the form (12), although not orthogonal in some Hilbert space, allow one to recover a bigger class of functions than in Point (a). The vector space V is not equal to L 2 (D, Leb), but the following result shows that it is much bigger than L 2 (D, µ). Proof To prove this result, it is more convenient to consider the variables (z, w) as in (10) instead of (x, y). The vector space V then becomes the set for some g 1 , g 2 ∈ L 2 ([a, b], Leb) and We first introduce the following notion of trace: we say that a function respectively. Our first claim is that any g ∈ L 2 ([a, b]×[−1, 1], Leb) which admits traces g 1 and g 2 at w = 1 and w = −1 respectively, belongs to V , and these traces are exactly the functions g 1 and g 2 in (18). To see that, we only have to check that Second, we claim that any g ∈ H 1 ([a, b] × [−1, 1], Leb) admits traces at w = 1 and w = −1 as defined above. Assume first that g ∈ C 1 ([a, b] × [0, 1]). Then, using the Cauchy-Schwartz inequality, Since in addition g(·, 1) L 2 (Leb) ≤ 4 g H 1 (Leb) by classical trace results (cf. e.g. [1, p. 196]), the function g → (g(·, 1), g − g(·, 1)) extends by density to a linear operator ψ continuous from Since obviously ψ 1 (g) + ψ 2 (g) = g, the claim is proved and the proof of Proposition 2.6 is completed.
Proof of Theorem 2.5 An eigenvector of A of the form H d (w)ψ(z) satisfies the Neumann boundary condition in D iff ψ ′ (a) = ψ ′ (b) = 0. The problem (13) with this boundary condition is a regular Sturm-Liouville problem with the weight ν defined in (14). Therefore (cf. e.g. [ , ν) composed of solutions to (13) on (a, b) with Neumann boundary conditions. All the corresponding eigenvalues are real, simple and the corresponding spectrum has no accumulation point. Now, we claim that forms an orthonormal basis of L 2 (D, µ). The orthonormal property follows from the fact that, if ϕ(x, y) = H d ( x−y x+y )ψ To prove that F is a basis of L 2 (D, µ), assume that f ∈ L 2 (D, µ) satisfies By the Cauchy-Schwartz inequality, for all d ≥ 2, the function In addition, Hence we deduce from Proposition 2.2 that this function is 0 for almost every z ≥ 0. Hence f = 0.
Thus F is an orthonormal basis of L 2 (D, µ) composed of eigenvectors of A. It is then classical to deduce that A admits no other eigenvector in this space, in the sense of point (a).
For point (b), let us first prove that the decomposition , Leb) and f 3 ∈ L 2 (D, µ). We only need to prove that this equality for as ε → 0. Therefore, This implies that b a f 2 2 (z)dz = 0, i.e. f 2 = 0. Similarly, f 1 = 0 and thus f 3 = 0. Since , the result then follows from the decomposition of f 1 and f 2 (resp. f 3 ) in the orthonormal basis {ψ k (x + y)} d≥2, k≥1 of L 2 (D, µ) ). To motivate the calculations of the next section, let us finally observe that, for all ϕ ∈ C 2 (R 2 + ), and that T L = L T .

On a family of bivariate polynomials
The goal of this section is to prove the existence of a family of polynomials in R of two variables X and Y , satisfying the family of relations for an integer d ≥ 0.
Before stating the main result of the section, let us recall some notation. R[X] is the set of polynomials on R with a single variable X and R[X, Y ] the set of real polynomials with two variables X and Y . The degree deg(P ) of a polynomial P ∈ R[X, Y ] is defined as the maximal total degree of each monomial of P . We define where only finitely many of the a i,j are nonzero. The real number a i,j will be called the (i, j)-coefficient of P .
For any P ∈ R[X, Y ] and for any d ≥ 0, we denote by [P ] d the sum of all monomials of P of degree d: In particular, [P ] d is homogeneous of degree d and P = ∞ i=1 [P ] i . We denote by ∆ i the first-order symmetrized discrete derivative with respect to the i-th variable: and by ∆ 2 i the symmetrized second-order discrete derivative with respect to the i-th variable: Note that the superscript in the notation ∆ 2 i does not correspond to the composition of the operator ∆ i with itself.
Finally, we define the linear operators on R[X, Y ] Then, adding and substracting the equations (20a) and (20b), the system (20) is equivalent to T P = 0.
We are going to prove the following result For any d ∈ {0, 2, 3, . . .}, the system (20) has a one-dimensional vector space of solutions. All these solutions, except 0, are of degree d. For d = 0, this is the vector space of constants, spanned by P 0 ≔ 1. When d ≥ 2, we denote by P d the unique solution to (20) It can be checked that the first polynomials are Before proving this result, let us give some properties of the polynomials P d , proved after the proof of Theorem 3.1.

Proposition 3.2
The polynomials P d , d ≥ 2, defined in Theorem 3.1 satisfy the following properties: where δ ij is the Kronecker symbol and by convention i j = 0 if j < 0 or j > i.

Proposition 3.3
For all d ≥ 2, the polynomial P d is given by the following formula: where In order to prove Theorem 3.1, we need the following lemma.
Proof of Lemma 3.4 To prove (25), it would be enough to expand both sides of the equation. We prefer to give a proof based on differential calculations, because it is related to the method used in the rest of the proof. First, let I = {1, 3, 5, . . .} and J = {2, 4, 6, . . .} be the sets of odd and even integers, respectively. Using the fact that, for Q ∈ R[X].
where Q (n) denotes the n-th derivative of Q, one has Since for all p, q ∈ N one easily checks that Now, for all n ∈ N, p∈J, q∈I, p+q=2n+1 This completes the proof of Lemma 3.4 (a).
To prove (b), let H p be the subspace of P d composed of all homogeneous polynomials of degree d. Recall the definition (19) of the operator T on has degree d. With this notation, it can be checked that , and similarly for the second variable. This yields Using the relation 4XY /(X + Y ) 2 = 1 − W 2 , we finally obtain that P ∈ Ker( T ) ∩ H d if and only if h solves (6). By Proposition 2.2, for all d = 1, this equation has a unique (up to a multliplicative constant) polynomial solution, which has degree d.
or, equivalently, Proof of Theorem 3.1 Fix d ≥ 0. We claim that, as a linear operator on P d , L is diagonalizable and its spectrum Sp P d (L) = {0, 1, . . . , d}. To see this, fix λ ∈ Sp P d (L) and P an eigenvector for this eigenvalue, with degree p. Writing as in the proof of Lemma 3.4 P = [P ] p + . . . + [P ] 0 , the equation LP = λP is equivalent to the fact that, for 0 ≤ n ≤ p Now, for any Q ∈ H k , one has Therefore, (29) for n = p imposes λ = p, and for n = p − 1, which allows one to compute recursively [P ] p−2 , . . . , [P ] 0 given any [P ] p ∈ H p and [P ] p−1 = 0. Since dim(H p ) = p + 1, the eigenspace corresponding to the eigenvalue p of L has dimension p + 1.
Let p ∈ {2, . . . , d} and assume that there exists which has a one-dimensional vector space of solutions in H p . Once [P ] p is fixed and since we have [P ] p−1 = 0, (31) can be used recursively to compute Proof of Proposition 3.2 Recall from the proof of Lemma 3.4 that a nonzero solution P of (20) for d ≥ 2 satisfies and Point (a) then follows from Proposition 2.2 and the value of the (d −1, 1) Observe that any polynomial solution of (6) with d ≥ 2 is divisible by (X− 1)(X + 1). As a consequence of the previous construction, any polynomial P such that (20) holds satisfies that XY divides [P ] d . Note also that [P ] d−1 = 0, which implies (b) by (31). Moreover, (31) also implies by induction that XY divides [P ] d−2k for all 0 ≤ k ≤ d/2, which yields the first part of (c). By Then it can be easily checked that, for all k ≥ 1, Therefore, (d) easily follows from (31) by induction.
Since this holds for any n, m ≥ 0, the same is true for all Q such that Q(X, Y ) = Q(Y, X). Now, any polynomial of the form P (X, Therefore, the fact that [P ] i is divisible by X − Y for i < d follows form (31) by induction. This ends the proof of (c).
As a consequence of (c), Point (g) is a simple consequence of points (e) and (f) and of formula (20a) with and applying (20b) to both terms on the l.h.s. yields This means that, for all k ≥ 2 and 2 ≤ d ≤ k, the vector It is straightforward to check that the matrix A k is self-adjoint for the inner product ·, · µ , where µ i = 1/i(k − i), which implies that two right eigenvectors of A k corresponding to different eigenvalues are orthogonal w.r.t. this inner product. This yields (23) for d = d ′ .
Finally, fix 2 ≤ d ≤ k. Using (20a) and (20b), we have Applying this equality inductively, we deduce that for some constant C. Now, by point (a) and Proposition 2.2, we have as k → +∞. Thus C = 2 and the proof of (h) is completed.
Proof of Proposition 3.3 In [14], the authors construct a family of functions of two variables satisfying relations close to (20), which they use to study neutral, multitype population processes with non-zero mutation or immigration. These functions are expressed in terms of the Hahn polynomials, defined for fixed parameters α > −1, β > −1 and N ∈ N by for all integer d ≥ 0. Karlin and McGregor proved that the rational function Let us define Passing to the limit in (34) and (35) proves that ψ d satisfies (20). Since ψ d is a polynomial, Theorem 3.1 entails (24). It only remains to check that The Chu-Vandermonde formula (cf. e.g. [8]) implies that which gives the expression of C d .

Spectral decomposition of neutral two-dimensional Markov chains
In this section, we consider neutral extensions of the two-dimensional birth and death chains in Z 2 described in the introduction. In the sequel, this family will be called N2dMC, for neutral two-dimensional Markov chains.
A N2dMC (X t , Y t ) t∈Z + is constructed by specifying first the Markov dynamics of Z t = X t + Y t in Z + . Assume that its transition matrix is where m≥0 p n,m = 1 for all n ≥ 1 and the state 0 is absorbing (p 0,0 = 1). Then, the process (X t , Y t ) t∈Z + is constructed as follows: if there is a birth at time t (i.e. if Z t+1 > Z t ), the types of the new individuals are successively picked at random in the population; if there is a death at time t (i.e. if Z t+1 < Z t ), the types of the killed individuals are successively picked at random in the population; finally, if Z t+1 = Z t , then X t+1 = X t and Y t+1 = Y t .
For example, the transition probability from (i, j) to (i + k, j + l) for (i, j) ∈ Z 2 + and k, l ≥ 0, k + l ≥ 1 is After some algebra, one gets the following formulas for the transition probabilities: for all l ≥ 0 and k ≥ 0 such that l + k ≥ 1, the Markov chain (X n , Y n ) n≥0 has transitions from (i, j) ∈ Z 2 + to with the convention that i j = 0 if i < 0, j < 0 or j > i. In particular, once one component of the process is 0, it stays zero forever. We denote by Π ≔ (π (i,j),(k,l) ) (i,j),(k,l)∈Z 2 + the transition matrix of the Markov process (X, Y ).
The state space of the Markov chain Z will be denoted by S Z , and the state space of (X, Y ) by S. We are going to consider two cases: the case where Z has finite state space S Z ≔ {0, 1, . . . , N} for some N ≥ 0, and the case where Z has infinite state space S Z ≔ Z + . In the first case, the state space of (X, Y ) is the set In the second case, S ≔ Z 2 + . We also define the sets S * ≔ S ∩ N 2 , T * N ≔ T N ∩ N 2 and S * Z ≔ S Z \ {0}. Finally, let Π 0 be the restriction of the matrix Π 0 to S * Z (i.e. the matrix obtained from Π 0 by suppressing the row and column of index 0) and let Π be the restriction of the matrix Π to S * .
Extending the usual definition for Markov chains, we say that a matrix M = (m ij ) i,j∈A is reversible with respect to the measure (µ i ) i∈A if µ i > 0 and µ i m ij = µ j m ji for all i, j ∈ A.
For all d ≥ 0, we define where we recall that the polynomials P d are defined in Theorem 3.1 and P 1 = P 1 . Note that, by Proposition 3.2 (e), a vector v i,j = P d (i, j)u i+j ∈ V d is characterized by the values of u k for k ≥ d only.
For all d ≥ 0, we also define the matrix Π d ≔ (p (38) All these notation, as well as those introduced in the rest of the paper, are gathered for convenience in Appendix A.
The following result is the basis of all results in this section.
for all d ≥ 0 and n ≥ 1, which can be written as for (X, Y ) ∈ Z 2 + . Similarly, an inductive use of (20b) yields Proposition 4.1 then easily follows from these equations and from the transition probabilities (36).

Eigenvectors of Π for finite state spaces
In the case where Z has finite state space, the main result of this section is the following.
is a right eigenvector for the eigenvalue θ of the matrix Π, where the polynomials P d are defined in Theorem 3.1 and where P 1 = P 1 . In addition, if d ≥ 2, v is also a right eigenvector for the eigenvalue θ of the matrix Π. (c) Assume that Π 0 admits a positive reversible measure (µ n ) n∈S * Z . Then, the matrix Π is reversible w.r.t. the measure and hence is diagonalizable in a basis orthonormal for this measure, composed of vectors of the form (39) for d ≥ 2.
In addition, Π is diagonalizable in a basis of eigenvectors of the form (39) for d ≥ 0.
Hence, the right eigenvectors of the transition matrix of a finite N2dMC can be decomposed as the product of two terms, one depending on each population size, but "universal" in the sense that it does not depend on the transitions matrix Π 0 of Z, and the other depending on the matrix Π 0 , but depending only on the total population size.

Remark 4.3
There is some redundancy among the right eigenvectors of Π of the form P (i, j)u i+j for P = P 0 , P (1) 1 or P 1 : if u is a right eigenvector of Π 1 , the vectors (iu i+j ) (i,j)∈S and (ju i+j ) (i,j)∈S are right eigenvectors of Π for the same eigenvalue. In particular, iu i+j + ju i+j is an eigenvector of Π of the form P 0 (i, j)u ′ i+j . This will also be true when Z has infinite state space.

Remark 4.4
In the following proof (and also in the case of infinite state space), no specific use is made of the fact that the matrix Π is stochastic. Therefore, Theorem 4.2 also holds true in the case of a continuous-time N2dMC, where the matrix Π 0 is now the infinitesimal generator of the process is composed of (N + 1) + N + (N − 1) + . . . + 1 = (N + 1)(N + 2)/2 = |T N | elements. Moreover, one can prove that it is linearly independent as follows: seen as a subspace of C N +1 by putting the d first coordinates to be zero. Given k ≤ N, the equality (41) for i + j = k combined with Therefore, F is a basis of C T N and, by point (a), the matrix Π has a Jordan normal form in this basis. Point (b) is then staightforward.
If Π 0 admits a positive reversible measure µ, it is straightforward to check that the vector ν in (40) is a reversible measure for Π, and hence the first part of Point (c) is true.
In addition, the matrix Π 1 is reversible w.r.t. the measure µ (1) n ≔ 2 n 2 µ n , n ∈ S * Z , which implies that Π 1 admits a basis of right eigenvectors orthonormal w.r.t. µ (1) . Similarly, Π 0 admits a basis of right eigenvectors orthonormal w.r.t. µ. By Point (a), this gives N +1 (resp. N) right eigenvectors of Π of the form (39) for d = 0 (resp. d = 1). Together with the basis of right eigenvectors of Π obtained above (extended by zero on {0} × N and N × {0}), this gives a basis of eigenvectors of Π and ends the proof of (c).
One may however recover a much larger class of matrices for which V d are stable vector spaces, by considering the algebra of matrices spanned by the matrices Π (n)± and Π (n) . Below, we study in detail such an example.
Consider an urn with N balls of three different colors and consider the following process: one picks a ball at random in the urn, notes its color, puts it back in the urn, picks another ball in the urn and replaces it by a ball of the same color as the first one. The number of balls of each colors then forms a Markov chain, which can be viewed as the embedded Markov chain of the 3-types Moran model, defined as follows: consider a population of N individuals, with 3 different types. For each pair of individuals, at rate 1, the second individual is replaced by an individual of the same type as the first individual in the pair.
Let i denote the number of balls of the first color, and j the number of balls of the second color. Then, there are N − i − j balls of the third color. The transition probabilities of this Markov chain are as follows: given that the current state of the process is (i, j), the state at the next time step is These transition probabilities do not not have the form (36). However, a variant of Proposition 4.1, and hence of Theorem 4.2, apply to this process, because of the following observation: let us construct the matrices Π + , Π − , Π + , Π + from the matrices Π + 0 = (p + n,m ) n,m∈S Z , Π − 0 = (p − n,m ) n,m∈S Z , Π + 0 = ( p + n,m ) n,m∈S Z , Π + 0 = ( p + n,m ) n,m∈S Z respectively, exactly as Π was constructed from Π 0 in (36), where , p + n,m = 0 otherwise. Then the transition matrix of the 3-colors urn model is given by In particular, the vector spaces V d for 0 ≤ d ≤ N are all stable for this matrix. The admits the reversible measure .
Hence the matrix Π admits a family of right eigenvectors null on the absorbing sets, which forms an orthonormal basis of L 2 (S * * , ν).
One easily checks, using (32), that v (i,j) = P d (i, j)u i+j is a right eigenvector of Π for the eigenvalue θ if and only if for all d ≤ k ≤ N where Now, the Hahn polynomials Q n (x; α, β, N) introduced in (33) satisfy (cf. e.g. [14]) Hence, for all 0 ≤ n ≤ N − d, (42) admits the polynomial (in k) solution of degree n If n ≥ 1, this polynomial must be divisible by (N − k), so we can define the polynomial R and the eigenvalue corresponding to the eigenvector P d (i, j) if n = 0, or Similarly as in the proof of Proposition 3.3, this family of eigenvectors can be seen as a singular limit case of those obtained in [14] for the multitype Moran model with mutation or immigration. Note that in the case of the 2-colors urn model, one can easily check that a basis of right eigenvectors of the corresponding transition matrix is given by (1, . . . , 1) and (NR Hence the spectrum is the same in the 2-and 3-colors urn models, although the multiplicity of each eigenvalue is different. In the case of two colors, the eigenvalues have the form 1 − k(k − 1)/N 2 for 0 ≤ k ≤ N, each with multiplicity 1 (except for the eigenvalue 1, with multiplicity 2). In the case of three colors, the eigenvalue 1 − k(k − 1)/N 2 has multiplicity k + 1 (except for the eigenvalue 1, which has multiplicity 3).
Concerning the eigenvectors in L 2 (S * * , ν). they are given by and for all 3 ≤ k ≤ N, the eigenspace for the eigenvalue 1 − k(k − 1)/N 2 is We shall end the study of this example by giving an apparently non-trivial relation between the polynomials P d and R (N,d) n . Because of the symmetry of the colors, we have

The case of infinite state space
The goal of this subsection is to extend Theorem 4.2 to the case where Z has infinite state space and the matrix Π 0 is a compact operator and admits a reversible measure. To this aim, we need first some approximation properties of Π by finite rank operators.

Approximation properties
Recall that Π 0 is a Markov kernel on Z + absorbed at 0 and that Π 0 denotes its restriction to N (equivalently, Π 0 is the sub-Markovian kernel on N corresponding to Π 0 with a Dirichlet condition at 0). We assume that Π 0 is reversible with respect to a positive measure µ on N (not necessarily finite). For any N ∈ N, consider the sub-Markovian kernel Π In other words, Π remains reversible with respect to µ, and thus, like Π 0 , it can be extended into a self-adjoint operator on L 2 (N, µ). We denote by ||| · ||| 0 the natural operator norm on the set of bounded operators on L 2 (N, µ), namely if K is such an operator, If furthermore K is self-adjoint, we have, via spectral calculus, The next result gives a simple compactness criterion for Π 0 . has finite range, if the above convergence holds, Π 0 can be strongly approximated by finite range operators and it is well-known that this implies that Π 0 is compact.
The converse implication can be proved adapting a standard argument for compact operators: assume that Π 0 is compact and let ε > 0 be fixed.
for some n < +∞ and ψ (1) , . . . , ψ (n) ∈ Π 0 (B), where B(ψ, ε) is the closed ball centered at ψ with radius ε. For any i ≤ n, since In other words, for all ϕ ∈ B, there exists i ≤ n such that Π 0 ϕ−Pr N i ψ (i) µ ≤ 2ε. This implies that In order to complete the proof, it only remains to check that If this was false, one could find a sequence (ϕ (N ) ) N ≥1 in B such that ϕ 0 . More generally, we consider two Markov kernels Π 0 and Π ′ 0 on Z + absorbed at 0, whose restrictions to N, Π 0 and Π ′ 0 , are both reversible with respect to µ. We associate to them Π and Π ′ defined on Z 2 + as in (36), and their respective restriction to N 2 , Π and Π ′ . We also define the matrices Π d and Π ′ d for d ≥ 1, as in (38). Note that Π and Π ′ are reversible with respect to ν, defined in (40) and it is straightforward to check that, for We will denote |||·||| and |||·||| d the operator norms in L 2 (N 2 , ν) and the sense that all the entries of this infinite matrix are non-negative), then and In particular, . We denote by v(d, u) the sequence v defined by the above r.h.s. The definitions of ν, µ (d) and Proposition 3.2 (h) enable us to see that the mapping is an isometry.
The Cauchy-Schwartz inequality and Proposition 3.2 (h) imply that v (d) ∈ L 2 (N d , µ (d) ) and since v is orthogonal to V ′ d , by the definition of V ′ d , the vector v (d) is orthogonal to L 2 (N d , µ (d) ), i.e. v (d) = 0. Fixing l ≥ 2 and applying Proposition 3.2 (f), one deduces from the equations v (d) l = 0 for 2 ≤ d ≤ l that v i,l−i = 0 for 1 ≤ i ≤ l − 1, and thus v = 0, ending the proof of (46). Now, Proposition 4.1 show that V ′ d is stable by Π and Π ′ and, more precisely,

It then follows from (46) that
Let us now come to the proof of (44). Let d ≥ 1 be fixed and define i.e.
For any u ∈ L 2 (N d+1 , µ (d+1) ), we get i |u i (x)|, and u i = 0 if i = d. It is clear that u ∈ L 2 (N d , µ (d) ), so, taking the supremum over u ∈ L 2 (N d+1 , µ (d+1) ), we have Hence, it only remains to show that ending the proof of (44). A similar computation using the fact that Again, let Π (N ) be the restriction of Π (N ) to T * N , which can thus be seen as a finite T * N ×T * N matrix. Similarly to the remark after the proof of Lemma 4.5, the spectrum of Π (N ) is the spectrum of Π (N ) plus the eigenvalue 0.

Spectral decomposition of infinite, compact, reversible N2dMC
The following result is an immediate consequence of Lemma 4.5 and Proposition 4.6.
Corollary 4.7 If Π 0 is compact and reversible, the same is true for Π.
We can now extend Theorem 4.2 to the infinite compact, reversible case.
which is a basis of the vector space in the sense that, for all v ∈ V , there exist unique sequences {α l } l≥1 , {β l } l≥1 , and {γ dl } d≥2, l≥1 and a unique a ∈ R such that where the series l α l ku (1),l k and l β l ku Example 1 Assume that Z n = X n +Y n is a birth and death process, i.e. that the matrix Π 0 is tridiagonal. Assume moreover that all the entries just above or below the diagonal are positive (except of course for p 0,1 , which is 0 since p 0,0 = 1). It is well-known in this case that there always exists a reversible measure µ for Π 0 . A well-known sufficient condition for the compactness of Π 0 is the case where this operator is Hilbert-Schmidt, which translates in our reversible, discrete case as i,j∈N For a birth and death process, letting p k (resp. q k ) denote the birth (resp. death) probability in state k and r k ≔ 1 − p k − q k ≥ 0, this gives As a side remark, note that Π is not necessarily Hilbert-Schmidt when Π 0 is, as the condition (i,j),(k,l) π (i,j),(k,l) π (k,l),(i,j) < ∞ is not equivalent to (50).
Proof Point (a) can be proved exactly the same way as Theorem 4.2 (a). The fact that compact selfadjoint operators admit an orthonormal basis of eigenvectors is classical. To prove (b), we only have to check that all these eigenvectors can be taken of the form (39) for d ≥ 2. This easily follows from the fact that Π d is compact selfadjoint in L 2 (N d , µ (d) ) for all d ≥ 2, from (46) and from Point (a).
The proof of (c) is similar to the proof of Theorem 2.5 (b). Fix v ∈ V and define a, v (1) , v (2) and v (3) as in (48). Since µ where the convergence of the first series holds for · µ (1) , and thus of the second series for · µ . A similar decomposition for v (2) and the use of (b) for v (3) complete the proof of (49). It only remains to observe that, for all v ∈ V , the equation i,0 = 0, one must have a = v 0,0 , v

On Dirichlet eigenvalues in the complement of triangular subdomains
In this section, we consider the same model as in the previous section, and we assume either that S = T N is finite or that S = Z 2 + and the restriction Π 0 of Π 0 to N is compact and reversible w.r.t. some measure µ. We recall that for all k ≥ 2. Note that S * 2 = S * . We also define S * 1 ≔ S \({0}×Z + ). Finally, for k ≥ 0, we call Π k the restriction of the matrix Π 0 to {i ∈ S Z : i > k}, and for k ≥ 1, Π k the restriction of the matrix Π to S * k . Note that this notation is consistent with the previous definition of Π 0 and that Π 2 = Π. Again, all the notations of this section are gathered for reference in Appendix A.

The case of finite state space
Let us first assume that S = T N for some N ≥ 1.
For 1 ≤ k ≤ N, the Dirichlet eigenvalue problem for the matrix Π in the set S * k consists in finding θ ∈ C and v in R T N , such that This is equivalent to finding a right eigenvector of Π k and extending this vector by 0 to indices in S \ S * k . For all k ≥ 1, we define θ D k as the supremum of the moduli of all Dirichlet eigenvalues in S * k . By Perron-Fröbenius' theory, θ D k is a Dirichlet eigenvalue in S * k . For all d ≥ 0, we also define θ (d) as the supremum of the moduli of all eigenvalues of Π corresponding to right eigenvectors of the form P d (i, j)u i+j . Again, by Theorem 4.2 (a) and Perron-Fröbenius' theory, θ (d) is an eigenvector of Π of the form P d (i, j)u i+j . Note that, for all d ≥ 2, because of Proposition 3.2 (c) and (e), any right eigenvector of Π of the form . The next result gives other inequalities concerning these eigenvealues. (a) Then, and If Π 0 is irreducible and p i,0 > 0 for some 1 ≤ i ≤ N, then θ (1) < 1.
Proof Since the matrix Π 0 is stochastic, it is clear from Theorem 4.2 (a) that θ (0) = 1. By (38), p n,m = m n p n,m for all 1 ≤ n, m ≤ N and thus Π 0 u = θu if and only if Π 1 v = θv, where v n = u n /n. Hence, the largest eigenvalue of . Ordering conveniently the states in S * 1 , the matrix Π 1 has the form for some rectangular nonnegative matrix Q, since the set {0, 1, . . . , N}×{0} is absorbing for the Markov chain. Again, since this matrix is block triangular, Since in addition Π N = p N,N Id and Π N = p N,N , Theorem 5.1 (a) will hold true if we prove that the sequences (θ (k) ) 1≤k≤N and (θ D k ) 2≤k≤N are both nonincreasing.
By Perron-Fröbenius' characterization of the spectral radius of nonnegative matrices (cf. e.g. [11]), for all 1 ≤ k ≤ N, where, by convention, the fraction in the r.h.s. is +∞ if u i = 0. Using the notation of the proof of Proposition 4.6, for all 1 ≤ k ≤ N − 1 and Taking the supremum over u ∈ R For all k ≥ 2, the Dirichlet eigenvectors in S * k belong to the vector space By Perron-Fröbenius' theory again, for all k ≥ 2, This ends the proof of (a).
In the case where Π k−1 is irreducible for some 2 ≤ k ≤ N − 1, it is clear that Π k and Π k are both irreducible. Then, by Perron-Fröbenius' theory, θ D k (resp. θ (k) ) is an eigenvalue of Π k (resp. Π k ) with multiplicity one, and the corresponding nonnegative eigenvector has all its entries positive. In addition, θ D k (resp. θ (k) ) is the only eigenvalue of Π k (resp. Π k ) corresponding to a positive eigenvector. In particular, the supremum in (58) (resp. (56) ) is attained only at vectors w ∈ V (k−1) having all coordinates corresponding to states in S * k−1 positive (resp. at vectors u ∈ (0, ∞) N k ). Hence, the inequalities in (59) and (57) are strict.
In the case where Π k−1 is irreducible for 3 ≤ k ≤ N − 1, let (u i ) i∈N k be a positive right eigenvector of Π k . Then the vector P k (i, j)u i+j belongs to U k and its restriction to S * k is a right eigenvector of Π k . However, its has positive and negative coordinates by Proposition 3.2 (g). Therefore (54) is proved.
Finally, if Π 0 is irreducible and p i,0 > 0 for some 1 ≤ i ≤ N, then the absorbing state 0 is accessible by the Markov chain from any initial state (possibly after several steps). It is then standard to prove that there exists n such that the sums of the entries of each line of ( Π 0 ) n is strictly less than 1. This proves that ( Π 0 ) n cannot have 1 as eigenvalue, and thus θ (1) < 1.

The case of infinite state space
Our goal here is to extend the previous result to the case where S = Z 2 + and Π 0 is compact reversible. So let us assume that Π 0 is compact and is reversible w.r.t. some measure µ.
In the case of infinite state space, when k ≥ 2 the Dirichlet eigenvalue problem for Π in S * k consists in finding θ ∈ C and v ∈ L 2 (Z 2 + , ν) (where the measure ν is extended by convention by zero on Z + × {0} ∪ {0} × Z + ) satisfying (52). Defining the vector space where Dirichlet eigenvectors are to be found k of the moduli of all Dirichlet eigenvalues in S * k is given by In view of Theorem 4.8 (c), the natural space to define the Dirichlet eigenvalue problem in S * 1 is where v (1) ∈ L 2 (N, µ) and v (3) ∈ L 2 (N 2 , ν) , (1) and v (3) are uniquely defined from v ∈ U 1 ). Then, the Dirichlet eigenvalue problem in S * 1 consists in finding θ ∈ C and v ∈ U 1 satisfying (52). We also define θ D 1 as the supremum of the moduli of all Dirichlet eigenvalues in S * 1 . For all d ≥ 1, we also define θ (d) as the supremum of the moduli of all eigenvalues of Π corresponding to right eigenvectors of the form P d (i, j)u i+j with u ∈ L 2 (N d , µ (d) ). By Theorem 4.8 (a), , ∀d ≥ 0.
In addition, using the notation V ′ d defined in the proof of Proposition 4.6, it follows from (46) that Comparing this with (60), we again deduce from Proposition 3.2 (c) and (e) that θ D d ′ ≥ θ (d) for all 2 ≤ d ′ ≤ d. Finally, since the matrix Π 0 is not reversible, we need to define θ (0) in a slightly different way: θ (0) is the supremum of the moduli of all eigenvalues of Π corresponding to right eigenvectors of the form P d (i, j)u i+j with u = a1 + v, where a ∈ R, 1 is the vector of R Z + with all coordinates equal to 1, v ∈ L 2 (N, µ) with the convention v 0 = 0.
Theorem 5.2 Assume that S = Z 2 + and that Π 0 is compact and reversible w.r.t. a positive measure µ. Then, for all d ≥ 1, θ D d is a Dirichlet eigenvalue of Π in the set S * d and θ (d) is a right eigenvalue of Π for an eigenvector of the form P d (i, j)u i+j with u ∈ L 2 (N d , µ (d) ). In addition, Proof For all k ≥ 1, the fact that θ (k) is an eigenvalue of Π for a right eigenvector of the form P k (i, j)u i+j with u ∈ L 2 (N k , µ (k) ) follows from Theorem 4.8 (a) and from the fact that Π k is compact reversible. Indeed, |||Π k ||| k is an eigenvalue of the compact reversible matrix Π k (the corresponding eigenvector can be obtained for example as the limit of nonnegative eigenvectors of finite dimensional projections of Π k ). The result follows since, necessarily, |θ (k) | ≤ |||Π k ||| k . For all k ≥ 2, we define Π (k) as the matrix whose restriction to S * k is Π k and with all other coordinates equal to zero. As the projection of a compact operator, this matrix is compact. Since it is trivially reversible for the measure ν, Theorem 4.8 (b) applies to Π (k) . Then, Proposition 4.6 tells us that θ D k = |||Π 2 is defined from Π (k) as Π 2 has been defined from Π. Therefore, the fact that θ D k is a Dirichlet eigenvalue for Π in S * k can be deduced exactly as above.
Recall the definition of Π (N ) in Section 4.2.1. For any N ∈ N, replacing Π by Π (N ) , we define similarly as above the quantities θ D,(N ) k and θ (k,N ) . Due to the remark after the proof of Proposition 4.6, we are again brought back to the finite framework. The following result is an immediate consequence of Lemma 4.5 and Proposition 4.6.
From this lemma and Theorem 5.1 immediately follow all the inequalities in Theorem 5.2 that concern θ D k for k ≥ 2 and θ (k) for k ≥ 1. As in the finite case, we easily deduce from the facts that p (1) i,j = j i p i,j and µ (1) ). Therefore, Since for all a ∈ R and v ∈ L 2 (N, µ), This yields θ (0) = 1 ≥ θ (1) .
In order to complete the proof, it only remains to check that θ D 1 = θ (1) and that θ D 1 is a Dirichlet eigenvalue in S * 1 . As in the finite case, Π 1 has the block-triangular form (55). Therefore, we obviously have θ D 1 ≥ θ D 2 . In addition, any Dirichlet eigenvalue in S * 1 which corresponds to an eigenvector in U 1 which is nonzero on the set of indices N × {0}, must be an eigenvalue of Π 0 corresponding to a right eigenvector in L 2 (N, µ). Now, if u ∈ L 2 (N, µ) satisfies Π 0 u = θ (1) u, then Π 1 v = θ (1) v with v i = u i /i and it follows from Theorem 4.8 (a) that the vector i i+j u i+j is a right eigenvector of Π 1 . Since this vector obviously belongs to U 1 , we obtain that θ D 1 = θ (1) and that θ D 1 is a Dirichlet eigenvalue in S * 1 .
6 Application to quasi-stationarity in nearly neutral finite absorbed two-dimensional Markov chains In this section, we restrict ourselves to the finite state space case for simplicity: let S = T N for some N ∈ N. The first coordinate will be referred to as type 1 and the second coordinate as type 2. Recall that the sets Z + × {0}, {0} × Z + and {(0, 0)} are absorbing for the N2dMC considered above, which means that each sub-population in the model can go extinct. This means that the transition matrix Π has the form where Q i (1 ≤ i ≤ 3) are square matrices and R i (1 ≤ i ≤ 2) rectangular matrices.
In this section, we study the problem of quasi-stationary distributions (QSD) and quasi-limiting distributions (QLD, see the introduction) for Markov processes, not necessarily neutral, whose transition matrix has the form (61-62). The classical case [6] for such a study is the case when Q is irreducible, which does not hold here. A general result is proved in Subsection 6.1. Our results of Section 5 are then applied in Subsection 6.2 to study the quasilimiting distribution of nearly neutral two-dimensional Markov chains.

Yaglom limit for general absorbing two-dimensional Markov chains
Let (X n , Y n , n ≥ 0) be a Markov chain on S = T N , with transition matrix of the form (61-62). We do not assume that this process is neutral. We call such processes A2dMC for "absorbed two-dimensional Markov chains". Under the assumption that the matrices Q 1 , Q 2 and Q 3 are irreducible and aperiodic, Perron-Frobenius' theory ensures the existence of a unique eigenvalue θ i with maximal modulus for Q i , which is real and has multiplicity one. Moreover, Q i admits unique (strictly) positive right and left eigenvectors for θ i , u i and v i respectively, normalised as v i 1 = 1 and v i u i = 1, where we use the convention that u i is a row vector and v i is a column vector, and where 1 denotes the column vector with adequate number of entries, all equal to 1. In the following result, we use the classical identification of column vectors and measures: for With this notation, using the results of Darroch and Seneta [6], v 1 ⊗ δ 0 and δ 0 ⊗ v 2 are trivial QSDs for the Markov chain (X, Y ). Theorem 6.1 Assume that the matrices Q 1 , Q 2 and Q 3 are irreducible and aperiodic, R 1 = 0 and R 2 = 0. Then, for any i ≥ 1 such that (i, 0) ∈ S, and similarly for the initial state (0, i) ∈ S, where L (i,j) denotes the law of the Markov chain (X n , Y n ) n≥0 with initial consition (i, j). Moreover, for any (i, j) ∈ S * , where and To give an interpretation of this result, we will say that there is extinction of type i conditionally on non-extinction if the QLD (64) gives mass 0 to all states with positive i-th coordinate. Conversely, we say that there is coexistence conditionally on non-extinction if the QLD (64) gives positive mass to the set S * . We also say that type 1 is stronger than type 2 (or type 2 is weaker than type 1) if θ 1 > θ 2 , and conversely if θ 2 > θ 1 .
Theorem 6.1 says that the limit behaviour of the population conditionally on non-extinction essentially depends on whether the largest eigenvalue of Q is θ 1 , θ 2 or θ 3 . If either θ 1 > θ 2 and θ 1 ≥ θ 3 or θ 2 > θ 1 and θ 2 ≥ θ 3 , the QLD is the same as if there were no individual of the weaker type in the initial population, and there is extinction of the weaker type conditionally on nonextinction. If θ 3 > θ 1 , θ 2 , there is coexistence of both types conditionally on non-extinction. Finally, when θ 1 = θ 2 , both types can survive under the QLD, so none of the types go extinct (according to the previous terminology), but there is no coexistence, as one (random) type eventually goes extinct. Observe also that the case θ 1 = θ 2 > θ 3 is the only one where the QLD depends on the initial condition.
Note also that, in the case where θ 3 ≤ max{θ 1 , θ 2 } and θ 1 = θ 2 , the QLD does not depend on any further information about the matrix Q 3 . In other words, if one knows a priori that θ 3 ≤ max{θ 1 , θ 2 } and θ 1 = θ 2 , the precise transition probabilities of the Markov chain from any state in S * have no influence on the QLD. The QLD is only determined by the monotype chains of types 1 and 2.
Our next result says that, for any values of θ 1 , θ 2 and θ 3 , all the QSDs of the Markov chain are those given in the r.h.s. of (64), when they exist and are nonnegative. Proposition 6.2 Under the same assumptions and notation as in Theorem 6.1, the set of QSDs of the Markov chain is composed of the probability in the case where θ 3 > max{θ 1 , θ 2 }.
In the case where θ 3 ≤ θ i for i = 1 or 2, we deduce from (69) that (θ 3 − θ i )µ i u i = av 3 R i u i . This is impossible since the l.h.s. of this formula is non-positive and the r.h.s. is positive as R i = 0, v 3 > 0 and u i > 0.
Similarly, it follows from (73) and Lemma 6.3 that The proof is easily completed in each of these cases.
Proof of Lemma 6.3 Assume first that θ 1 > θ 3 . One easily checks that

The nearly neutral case
Since Π in (61) is a block triangular matrix, we have where Sp ′ (A) denotes the spectrum of the matrix A, where eigenvalues are counted with their multiplicity.
In the case of a N2dMC satisfying the assumptions of Theorem 6.1, with the notation of Section 4, we have Q 1 = Q 2 = Π 0 . By Theorem 4.2 (b) and Remark 4.3, {1}∪Sp ′ (Q 1 )∪Sp ′ (Q 2 ) is the set of eigenvalues corresponding to right eigenvectors of Π of the form P 0 (i, j)u i+j and P 1 (i, j)u i+j , counted with their multiplicity. More precisely, Sp ′ (Q 1 ) corresponds to eigenvectors of the form P (1) 1 (i, j)u i+j , and Sp ′ (Q 2 ) to eigenvectors of the form P 1 (i, j)u i+j . In particular, θ 1 = θ 2 = θ (1) = θ D 1 , with the notation of Theorem 5.1. Moreover, since Q 3 = Π (2) , Theorem 5.1 shows that θ 3 = θ (2) = θ D 2 < θ 1 = θ 2 and Sp ′ (Q 3 ) is the set of eigenvalues corresponding to right eigenvectors of Π of the form P d (i, j)u i+j for d ≥ 2, counted with their multiplicity.
In other words, with the terminology defined after Theorem 6.1, coexistence is impossible in the neutral case. Since the eigenvalues of Q 1 , Q 2 and Q 3 depend continuously on the entries of these matrices, we deduce that coexistence is impossible in the neighborhood of neutrality: Corollary 6.4 Let Π be the transition matrix of some fixed N2dMC in T N such that Π 0 and Π 1 are both irreducible and there exists i ∈ {1, . . . , N} such that p i,0 > 0. For any A2dMC (X, Y ) in T N with transition matrix Π ′ sufficiently close to Π, coexistence is impossible in the QLD of (X, Y ). Let θ ′ 1 , θ ′ 2 , θ ′ 3 denote the eigenvalues θ 1 , θ 2 , θ 3 of Theorem 6.1 corresponding to the matrix Π ′ . If θ ′ 1 = θ ′ 2 , the QLD of (X, Y ) is the trivial QSD corresponding to the stronger type: if θ ′ 1 > θ ′ 2 , the QLD of (X, Y ) is v ′ 1 ⊗ δ 0 , where v ′ 1 is the QLD of (X, 0), and if θ ′ 2 > θ ′ 1 , the QLD of (X, Y ) is δ 0 ⊗ v ′ 2 , where v ′ 2 is the QLD of (0, Y ).

A Notations
We gather here all the notations used at several places in the paper. Most of these notations are introduced in Sections 4 and 5.

A.1 General definitions
• For any measurable subset Γ of R d and any σ-finite positive measure µ on Γ, L 2 (Γ, µ) is the set of Borel functions f on Γ defined up to µ-negligible set such that Γ f 2 dµ < +∞. We denote by ·, · µ the canonical inner product on L 2 (Γ, µ) and · µ the associated norm. In the case when Γ is discrete, we make the usual abuse of notation to identify the measure µ and the corresponding function on Γ.
• For all I × I square matrix M, where I is a finite or denumerable set of indices, and for all J ⊂ I, we call "restriction of the matrix M to J" the matrix obtained from M by removing all rows and columns corresponding to indices in I \ J.

A.4 Matrices
Π 0 = (p n,m ) n,m∈S Z is a stochastic matrix such that p 0,0 = 1. Π k is the restriction of Π 0 to the set of indices S Z ∩ N k+1 , for all k ≥ 0. Pr N is the projection operator on R N defined by Pr N (u 1 , u 2  is the Markovian kernel on Z + whose restriction to N is Π A.5 Measures and vectors µ = (µ i ) i∈S * Z is a reversible measure for the matrix Π 0 . ν = (ν (i,j) ) (i,j)∈S * , where ν (i,j) = (i + j)µ i+j ij .
A.6 Operator norms (in the infinite, reversible case) |||·||| 0 is the natural operator norm on the set of bounded operators on L 2 (N, µ). |||·||| is the natural operator norm on the set of bounded operators on L 2 (N 2 , ν). ||| · ||| d is the natural operator norm on the set of bounded operators on L 2 (N d , µ (d) ), for all d ≥ 1.

A.7 Eigenvalues
We refer to Sections 5.1 and 5.2 for precise definitions in the finite and infinite cases.
θ D k is the biggest Dirichlet eigenvalue of Π in S * k , for all k ≥ 1. θ (d) is the biggest eigenvalue of Π corresponding to right eigenvectors of the form P d (i, j)u i+j , for all d ≥ 0.

A.8 Vector spaces
Finite case U k = {v ∈ R T N : v = 0 on T N \ S * k } for all k ≥ 1.