Analyzing Generalized Pólya Urn Models using Martingales, with an Application to Viral Evolution

Consider the following generalized Pólya Urn model: An urn starts out with $u$ white balls and $v$ black balls, with $u+v>0$. At each step $i=1,2,3,\dots$, a ball in the urn is chosen uniformly at random. If the chosen ball is white, a black ball is added to the urn with probability $p_{W}$ and a white ball is added with probability $1-p_{W}$, where $0<p_{W}<1$. If the chosen ball is black, a white ball is added to the urn with probability $p_{B}$ and a white ball is added with probability $1-p_{B}$, where $0<p_{B}<1$. The originally chosen ball is then returned to the urn, so that at step $i$, the total number of balls equals $u+v+i$. This model, known as the randomized play-the-winner (RPW) model, has been widely studied in theoretical and applied contexts, with applications ranging from clinical trials to genetic mutations. Several papers, including Wei and Durham (1978), Smythe and Rosenberger (1995), Smythe (1996), and Rosenberger and Sriram (1997) have considered the RPW model’s asymptotic properties, deriving limit theorems related to the asymptotic fraction of white (or black) balls in the urn.

The distribution of this fraction after finitely many steps, however, is less well understood. Rosenberger and Sriram (1997) prove an expression for its expectation, which we re-derive using an alternative approach. Matthews and Rosenberger (1997) propose an expression for its variance; however, we claim that this expression is erroneous, as it returns the incorrect variance for the special case $p_{W}=p_{B}=\frac{1}{2}$, which reduces the RPW model to the Binomial$(\frac{1}{2})$ model (see Appendix). In this paper, we introduce a new approach to compute the variance of the number of white balls in the RPW model after finitely many steps. Our method involves rewriting the RPW model as a martingale process—a transformation that may be applied not only to the RPW model, but to any process represented by a sequence of random variables $M_{0},M_{1},M_{2},\dots$ with finite first and second moments such that ${\mathbb{E}}[M_{i+1}|M_{i}]$ is a non-constant linear function of $M_{i}$ and $\text{Var}[M_{i+1}|M_{i}]$ is a quadratic function of $M_{i}$, for all $i$. We then use our variance formula to approximate the full probability mass function (PMF) of $M_{i}$ when $p_{W}$ and $p_{B}$ are small. For such $p_{W}$, $p_{B}$, the RPW model aptly characterizes branching processes that arise in viral genetics, with the balls representing viral particles and colors representing variants. We conclude with a novel method of estimating mutation rates of pathogens by fitting our approximate PMF to genome sequencing data of viruses.

Our central idea behind computing the variance of any $M_{n}$ in the aforementioned process $\{M_{i}\}_{i\geq 0}$ is to construct a martingale $(X_{i},{\mathcal{F}}_{i})$ such that each $X_{i}$ is a linear function of $M_{i}$ (and ${\mathcal{F}}_{i}$ is the canonical $\sigma$-algebra, $\sigma(X_{0},\dots,X_{i})$). For a martingale, computations of ${\mathbb{E}}[X_{n}]$ and $\text{Var}[X_{n}]$ are more straightforward, and from them, we may obtain ${\mathbb{E}}[M_{n}]$ and $\text{Var}[M_{n}]$. Our first step in developing this method is an elementary lemma that allows us to compute ${\mathbb{E}}[M_{n}]$ for any $n$, when ${\mathbb{E}}[M_{i+1}|M_{i}]$ is a non-constant linear function of $M_{i}$: