Let's play a game. I start with $N$ indistinguishable tokens, and I wait $T$ turns. Every turn, each token has probability $p$ of disappearing. I want an analytic formula for the entropy of this process, as a function of $N$, $T$, and $p$.
The calculation is straightforward for $N=1$ and $T=\infty$. The probability $p_i$ that our (only) token disappears at turn $i$ is $(1-p)^{i}p$, and the entropy $E$ is given by:
$E(N=1, T=\infty, p=p) = \sum_{i}^{\infty}p_i\ln(p_i)$
$=\sum_{i}^{\infty}(1-p)^{i}p\ln((1-p)^{i}p)$
$=p\ln(1-p)\sum_{i}^{\infty}i(1-p)^{i} + p\ln(p)\sum_{i}^{\infty}(1-p)^{i}$
$=p\ln(1-p)\frac{1-p}{p^2} + p\ln(p)\frac{1}{p}$
$=\frac{1-p}{p}\ln(1-p) + \ln(p)$
For $N=2$ and $T=\infty$, my calculation (not shown) is a lot uglier, but simplifies down to:
$E(N=2, T=\infty, p=p) = \frac{2-2p}{2-p}\ln(2) + 2\ln(p) + \frac{2-2p}{p} \ln(1-p)$
I'm about to calculate the $N=3$, $T=\infty$ case, but I've got the feeling I'm reinventing the wheel. Is the formula for $E(N, T, p)$ known? I'm particularly interested in the $T=\infty$ case. A good approximation is almost as useful to me as an exact formula, but I'm interested in both small and large values of $N$.
As a sanity check, we can compare against a Python simulation:
#!/usr/bin/python3
import numpy as np
from scipy.stats import binom
def entropy(n_initial, p, n_steps, n_trials):
# Simulate many stochastic binomial decays, return the average of
# the log of their binomial "penalty"
log_penalty = np.zeros(n_trials, dtype=np.float64)
n = n_initial * np.ones(n_trials, dtype=np.int64)
for i in range(n_steps):
num_losses = np.random.binomial(n, p)
log_penalty += binom.logpmf(k=num_losses, n=n, p=p)
n -= num_losses
return log_penalty.mean()
p=0.1
print(entropy(n_initial=1, p=p, n_steps=1000, n_trials=10000))
print(((1-p)/p)*np.log(1-p) + np.log(p))
print(entropy(n_initial=2, p=p, n_steps=1000, n_trials=10000))
print((2-2*p)/(2-p)*np.log(2) + 2*np.log(p) + ((2-2*p)/p) * np.log(1-p))
Please forgive/correct me if I've made errors in my math or I'm using the wrong terms; I'm an experimental physicist, not a mathematician, and my formal math is rusty.
