Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

So I'm looking to calculate the probability that a bacterial population dies out with the following conditions:

  • Initial population is 1
  • At each iteration the bacteria can die, do nothing, divide into 2 or divide into three. So each bacteria can change into either 0, 1, 2, 3 bacteria.
  • Each event has a 25% likelihood.

So I know from simulation that the likelihood of the population dying out is ~36.3%. However, I'm looking for a method to calculate this analytically.

Some MATLAB code for the simulation:

ntries = 100000;
mnum = 10000;
res = ones(ntries, 1);

for i = 1:ntries
    while res(i) > 0 && res(i) < mnum
        res(i) = sum(randi(4,res(i), 1)-1);



share|cite|improve this question
Are you sure it's not 41.4%? – joriki Apr 28 '11 at 20:24
Nope I reliably get 0.3632 although this assumes that after a population gets to 10,000 then it won't die out ... so that might be where I'm loosing 5%. – JudoWill Apr 28 '11 at 21:25
@JudoWill: Three people have independently arrived at 41.4%, so I daresay it's likely that the mistake is on your side :-) – joriki Apr 28 '11 at 21:28
@JudoWill: I've never used MATLAB, but I think the "res(i)" after "sum" shouldn't be there -- you want to sum up the random values; there's no reason to multiply them with the current number of bacteria. – joriki Apr 28 '11 at 21:36
yup, my bad ... I get ~41% now – JudoWill Apr 28 '11 at 21:37
up vote 3 down vote accepted

Call $q_n$ the probability that the initial bacteria dies out at or before generation $n$ and condition by the size $i$ of the first generation. The event that every bacteria of this generation dies out at or before generation $n+1$ has probability $q_n^i$ and $i=0$, $1$, $2$ or $3$ with probability $1/4$, hence $q_{n+1}=g(q_n)$ with $$ g(q)=q^0/4+q^1/4+q^2/4+q^3/4. $$ Since $q_0=0$, $q_n=g\circ g\circ\cdots\circ g(0)$ ($g$ composed $n$ times). One can check, for instance by considering the graph of $g$ on the interval $(0,1)$, that $q_n$ converges to the smallest nonnegative root $q$ of the equation $q=g(q)$.

Now $g(q)=q$ if and only if $0=q^3+q^2-3q+1=(q-1)(q^2+2q-1)$, hence $q=1$ or $(q+1)^2=2$, that is $q=-1\pm\sqrt{2}$. Finally the population dies out with probability $q=\sqrt2-1\approx41\%$.

See this chapter of lecture notes by Steve Lalley, section 1.2.

share|cite|improve this answer
Either .41 or 41% :-) (Also $\approx$, not $=$) – joriki Apr 28 '11 at 20:26
@joriki: Yep. :-) Thanks. – Did Apr 28 '11 at 20:29

I find the probability to be $\sqrt{2}-1\approx 41.4\%$. If we let $p(n)$ be the probability that a population of $n$ dies out, we have $p(0)=1, p(n)=(p(n-1)+p(n)+p(n+1)+p(n+2))/4$. I first just did a numeric relaxation in Excel and recognized $p(1), p(2)=3-\sqrt{2}, p(3)=5\sqrt{2}-7.$ This was aided by the inverse symbolic calculator. The recurrence becomes $p(n+2)=-p(n+1)+3p(n)-p(n-1)$, whose characteristic polynomial has roots $1$, $\sqrt{2}-1$, and $-1-\sqrt{2}$. Finally we can check that $\sqrt{2}-1=(1+\sqrt{2}-1+3-2\sqrt{2}+5\sqrt{2}-7)/4$

share|cite|improve this answer

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.