On a species survival model

In this paper we provide some sharp asymptotic results for a stochastic model of species survival recently proposed by Guiol, Marchado, and Schinazi.

thus based on general properties of Markov chains, and suitable variations of the result can in principle be carried out to other similar models.
In this paper we focus on the model introduced in [7] as is. Our aim is to elucidate the underlying mechanism responsible for the phenomenon described in [7] by sharpening the estimates that lead to the major qualitative statement therein. We proceed with a description of the Guiol, Marchado, and Schinazi (GMS) model.
In contrast to the Bak-Sneppen model, the number of species in the GMS model is random and changes in time, and only the species with the lowest fitness is randomly replaced. The local interaction between species is not considered in the GMS model, and therefore the spatial structure of the population is of no importance. Let p > 1 2 be given and denote q = 1 − p. Let Z + denote the set of non-negative integers and let X = (X n : n ∈ Z + ) be a discrete-time birth and death process with the following transition probabilities: from each site, X n increases by 1 with probability p; from each site different than 0, X n decreases by 1 with probability q = (1 − p); finally, at 0, X n stays put with probability q. Thus X is a nearest-neighbor transient random walk on the integer lattice Z + with holding times and reflection barrier at zero. Throughout the paper we assume that X 0 = 0 with probability one.
The model describes the evolution of a population of species: newborn species are given "equal opportunity" by assigning them a random fitness, and the least fit species are the first to die. A jump to the right represents birth of a new species, whereas a jump to the left represents death of an existing species. Thus X n represents the number of species alive at time n. When a new species is born, it is assigned a fitness. The fitness is a uniform [0, 1] random variable independent on the fitness of all previously born species as well as of the path of the process X. When X jumps to the left, the species with the lowest fitness is eliminated. We remark that, in a different context, a similar model was considered by Liggett and Schinazi in [9].
Fix f ∈ (0, 1). We examine the model by considering two coupled random processes, L = (L n : n ∈ Z + ) (for lower or left) and R = (R n : n ∈ Z + ) (respectively, for right), where L n denotes the number of species alive at time n whose fitness is less than f while R n denotes the number of the remaining species alive at time n.
Observe that L n increases by 1 if X n does and the newborn species has fitness less than f, and L n decreases by 1 whenever X n decreases by one and L n is not zero. The value of L n remains unchanged when either X n increases by 1 and the newborn species has fitness at least f or X n decreases by 1 and L n = 0.
Notice that when it is not at zero, the process L evolves as a nearest-neighbor random walk with probability pf of jumping to the right, probability q of going to the left, and probability 1 − pf − q of staying put. When at zero, L n jumps to the right with probability pf , and stays put with the complementary probability 1 − pf. Thus L is a Markov chain. Since p > q, it is positive recurrent if pf < q, null-recurrent if pf = q, and is otherwise transient. In what follows we will denote the critical value q/p of the parameter f by f c .
Consider the process B = (B n : n ∈ Z + ), where B n is the total number of species born by time n with fitness at least f . Observe that B is a non-decreasing Markov chain (formed by sums of i.i.d. Bernulli variables), which jumps one step up with probabilities p(1 − f ) = p − q and stays put with the complementary probability 2q. It is shown in [7] that when f = f c , we have that for any ε > 0 while each species with fitness less than the critical value disappears after a finite (random) time. It follows from (1) that the distribution of species with fitness higher than f c approaches a uniform law. We sharpen the above result by proving the following theorem. Recall that f c = q/p.
where N (0, 1) denotes a mean-zero Gaussian random variable with variance one, and ⇒ stands for convergence in distribution.
3. (L n : n ∈ N) is a recurrent Markov chain (visiting zero infinitely often) and hence each species with fitness less than the critical value becomes eventually extinct with probability one.
The proof of the theorem is included in Section 2. Notice that by the law of large numbers B n ∼ n(1 − f c ), a.s., as n → ∞. Moreover, the fitness of any randomly chosen species from B n is distributed uniformly on (f c , 1). Hence the theorem implies that this uniform law on (f c , 1) is the limiting distribution for the fitness of a randomly chosen species being alive at time n and having fitness larger than f c .

Proof of Theorem 1
There is no loss of generality assuming that X is obtained recursively from an i.i.d. sequence of Bernoulli random variables J = (J n : n ∈ Z + ) with P (J n = 1) = q, P (J n = 0) = p, as follows: X 0 = 0 and for n ∈ Z + we have Let G = (G n : n ∈ Z + ) be a Markov chain formed by the pairs G n = (L n , J n ). Figure  1 illustrates the transition mechanism of G.
2.1 Reduction from B n − R n to an occupation time of G Let ∆ = (∆ n : n ∈ Z + ) be the process defined through ∆ n = B n − R n . Notice that ∆ n increases by 1 if and only if R n decreases by 1, and otherwise stays put. That is Thus by (2), ∆ n+1 − ∆ n = (1 − s n )J n 1 {Ln=0} , and we have Figure 1: Transition probabilities for G Observe that η n is the occupation time (number of visits) of G at (0, 1) up to time n−1 and that n−1 i=0 s i J i is bounded above by ∞ i=0 s i . Since X is transient and X 0 = 0, the latter series can be represented in the form where g is the total number of visits to zero (distributed geometrically, according to P [5, p. 274]) and h i are successive holding times of X at zero (i.i.d. sequence of geometric random variables, independent of g, with P (h i = k) = q k−1 p, k ∈ N). In particular, Since L and consequently G are recurrent (and thus ∆ n is a non-decreasing sequence converging to +∞ with probability one), it therefore suffices to show that Theorem 1 holds with ∆ n = B n − R n replaced by η n in its statement.
Excursion decomposition for the path of L. We decompose the path of L into a sequence of successive excursions from 0, each one lasting until (but not including) the first time L returns to 0 from 1. Observe that once the Markov chain L is at zero, it will stay put until J k = 0 and the fitness of the newborn particle is less than f c . Therefore, the holding time of L at zero during one excursion is a geometric random variable with parameter pf c = q. Leth k be the i.i.d. sequence of holding times at zero during successive excursions of L from zero. That is, P (h k = n) = p n−1 q, n ∈ N.

Define
Notice that (µ k : k ∈ N) is an i.i.d. sequence and We now compute µ := E(µ k ) using the fact that the sequence of pairs G n = (J n , L n ) forms a Markov chain, the transition mechanism of which is illustrated in Fig. 1. The value of µ is equal to the expected number of visits by this Markov chain to the state (0, 1) during the period of time starting at the state (0, 1) with probability q and at (0, 0) with probability p, and lasting until G leaves the set {(0, 0), (0, 1)}. We thus have, using first step analysis, Hence µ = 1. Consequently, using (3) and the law of large numbers, we obtain η n ∼ N n as n → ∞, a.s. (4) Reduction to a simple random walk The duration of the k-th excursion of L from zero is equal to the sumh k + α k of the holding time at zeroh k distributed as Geom(q), and the time α k that it takes for L to return from 1 to zero. Notice that For m ∈ N let V m = m k=1 (h k + α k ) be the total duration of the first m excursions of L from 0. Then N k = max{m : V m ≤ k}. With each excursion we can associate a skeleton, which is the path obtained from the excursion by omitting all transitions from a state to itself, for all states. If we let τ k denote the length of the skeleton, then due to the choice of f c it follows that τ k has the same distribution as the time required for the simple symmetric random walk on Z to get back to 0 starting from 0. Furthermore, the time spent by L at each visit to a site is a geometric random variable, Geom(2q) for sites different than 0 and Geom(q) for 0. Therefore, the length of a single excursion itself is a sum of one Geom(q) random variable plus a sum of τ k − 1 independent Geom(2q) random variables. If we replaceh k with a Geom(2q), then the resulting modified "excursion time" becomes a sum of τ k copies of a Geom(2q) random variable. Let V m denote the total length of the first m excursions modified in this way. Let N k denote the number of such excursions occurred by time k, that is N k = max{m : V m ≤ k}. Then V m stochastically dominates V m . By the law of large numbers, where h 1 is a geometric random variable with parameter 2q. Letting T m = m k=1 τ k , we obtain where (h k : k ∈ N) is an i.i.d. sequence of random variables, each one distributed as Geom(2q). Thus, by the law of large numbers, as m → ∞, a.s. (6) Notice that T m is distributed the same as the total length of the first m excursions from zero of a simple symmetric random walk.
2.2 Completion of the proof: CLT and LIL for η n LIL for η n We need the following result. Although the claim is a "folk fact", we give a short proof for the sake of completeness. Since φ −1 (x) ∼ x 2 /(2 ln ln x) as x → ∞, we obtain lim inf n→∞ γ n n 2 /(2 ln ln n) = 1.
Observe that γ n = n−1 i=0 Y i . The distribution of t 1 , the time to return to 0 starting from 0 for S, is equal to the distribution of t 1 , the time to return to 0 starting from 0 for the reflected random walk |S|. Since t 1 is equal in distribution to 1 + Y 1 , lim inf m→∞ T m m 2 /(2 ln ln m) = lim inf n→∞ γ n + n n 2 /(2 ln ln n) = 1, completing the proof of the lemma.
Using the lemma along with (5) and (6) CLT for η n We now turn to the proof of the central limit theorem. It is well known (see for instance [5, p. 394]) that lim m→∞ E e −θTm/m 2 = E e − √ 2θ , θ ≥ 0.
Therefore, it follows from (5) and (6)  The function θ → e −c √ 2θ , θ ≥ 0, is the Laplace transform of a positive stable law with index 1/2 whose density function is given by (see for instance [5, p. 395 We intend to use this formula with c = 1 √ 2q . Set V 0 = 0 and observe that for all k ∈ N, Fix s > 0 and let u = √ ks. Then Differentiating the right-hand side we obtain Therefore N k / √ k converges weakly to the absolute value of a centered normal random variable with variance equal to c −2 = 2q. Combining this with (4) completes the proof of the central limit theorem for η n .