Imagine there are 10 000 Thanosses on a planet, each having an infinity gauntlet. Each Thanos snaps its fingers ones, removing every other remaining Thanos with a change of 50%. After all Thanosses have either snapped there fingers or are disintegrated, what is the expected value of Thanosses remaining?
This riddle was inspired by the following riddle: https://fivethirtyeight.com/features/how-many-earthlings-would-survive-63-thanos-snaps/
You can watch the following video made with a horrible automatic voice, or read further to see the same argument layed out in text.
There are two kinds of Thanosses. One kind of Thanos is a Thanos that still has to snap it's fingers, in other words, that has a working gauntlet. Let us call this kind of Thanos an active Thanos. The other kind of Thanos is the Thanos that has already snapped it's fingers. It has no working gauntlet, and is therefore harmless. We call this kind of Thanos a passive Thanos.
If we forget for a moment the Thanos that is currently snapping, each other thanos has a chance of exactly one half to be disintegrated. This can be modelled by each Thanos flipping a coin.
- Head means Thanos survives
- Tail means Thanos dies
The distribution of the number of heads, and by consequence, the number of remaining Thanosses, can be described by the binomial distribution, with p being $1/2$.
Now, for each snap, there are 3 different kinds of Thanosses,
- Snapping Thanos (1 active): chance of 1 to have 1 passive Thanos
- Other active Thanosses ($a-1$ active): chance of $$\binom{k}{a-1}\frac{1}{2^{a-1}}$$ to have $k$ active Thanosses
- passive Thanosses ($p$ passive): chance of $$\binom{l}{p}\frac{1}{2^{p}}$$ to have $l$ active Thanosses
We can combine these equations to write what the chances are to remain with exactly k active and l plus one passive thanosses.
$$ \binom{k}{a-1}\frac{1}{2^{a-1}} \binom{l}{b}\frac{1}{2^{b}}$$
We can also look at what happens if we add an extra passive Thanos. This new passive Thanos does not affect the survival chance of any of the other Thanosses. It's survival chance is only affected by the number of remaining active Thanosses. This means that we can split the system in a part dependant on the active and on the passive Thanosses.
$$E(a,p+1) = E(a,p) + P(a)$$
where $P(a)$ is a single passive Thanos's survival chance.
Combining this with the fact that the expected value equals the sum of each of the chances times the expected value of each corresponding system gives us a nice recursive formula for the expected value.
$$E(a,p) = \sum_{k=0}^{a-1} \sum_{l=0}^{p} \binom{k}{a-1}\frac{1}{2^{a-1}} \binom{l}{p}\frac{1}{2^{p}} \left(E(k,l+1)\right)$$
where $E(0,p) = p$
Both formulas can be combined, giving us this nice expression:
$$A(a)+pP(A) = \sum_{k=0}^{a-1} \sum_{l=0}^{p} \binom{k}{a-1}\frac{1}{2^{a-1}} \frac{1}{2^{p}} \left(A(k)+(l+1)P(k)\right)$$
This can be solved to A by plugging in zero
$$A(a) = \sum_{k=0}^{a-1} \binom{k}{a-1}\frac{1}{2^{a-1}} \left(A(k)+P(k)\right)$$
We can also plug in one and substact A:
$$P(A)= \sum_{k=0}^{a-1} \binom{k}{a-1}\frac{1}{2^{a-1}} frac{1/2} P(k)$$
If we multiply the first equation with two and adding the second equation we get:
It can easily be seen that this expression remains constant, but we are going to proof it with string induction:
To proof: $A(a)+2P(a) = A(0)+2P(0) $
Base case: Trivial
Assume $A(k)+2P(k)$ is constant for $k \in [0,a-1]$:
$$A(a)+2P(a) = \sum_{k=0}^{a-1} \binom{k}{a-1}\frac{1}{2^{a-1}} \left(A(k)+2P(k)\right)$$
$$= \left(A(0)+2P(0)\right)\sum_{k=0}^{a-1} \binom{k}{a-1}\frac{1}{2^{a-1}} $$
$$= \left(A(0)+2P(0)\right) \times 1 $$
The value of the base case equals two, so the value of all other elements must also be two. This means that whenever there are two passive thanosses, the expected value will equal two. If there are less then two passive thanosses, like in our actual problem statement where there are zero, the number of remaining thanosses is less then two.
This simplifies our solution to the following equations.
$$P(a)= \sum_{k=0}^{a-1} \binom{k}{a-1}\frac{1}{2^a} P(k)$$
$$A(a)= 2-2P(a)$$
There is no nice formula for P, but it can be generated by a computer program. The following rust program generates all the values
use statrs::distribution::{Binomial, Discrete};
fn main() {
let mut p: Vec<f64> = vec![0.0; 100001 as usize];
p[0] = 1.0;
for n in 1usize..(1<<20) {
let binom = Binomial::new(0.5, (n-1) as u64).unwrap();
for k in 0usize..n {
p[n] += binom.pmf(k as u64) * p[k] * 0.5 ;
}
if vec![1,100,1000,10000,100000].iter().any(|x| *x == n) {
println!("{} & {} &{} \\\\",n,p[n],2.0-2.0*p[n]);
}
}
}
Running our program and taking the 10 000 result nicely gives us our solution of 1.999711521253186 which is very close to two, as expected.