The minimal observable clade size of exchangeable coalescents

TL;DR

This paper investigates the asymptotic behavior and moments of the minimal observable clade size in exchangeable coalescent models with mutations, providing insights relevant for genetic data analysis and model selection.

Contribution

It introduces asymptotic results and recursive formulas for the moments of the minimal observable clade size in $ ext{Lambda}$-coalescents with mutations, a quantity previously unobservable in real data.

Findings

Asymptotic behavior of $O_n$ as $n o $ derived.

Recursive formulas for all moments of $O_n$ established.

$O_n$ provides an upper bound for the minimal clade size.

Abstract

For $\Lambda$-$n$-coalescents with mutation, we analyse the size $O_n$ of the partition block of $i\in\{1,\ldots,n\}$ at the time where the first mutation appears on the tree that affects $i$ and is shared with any other $j\in\{1,\ldots,n\}$. We provide asymptotics of $O_n$ for $n\to\infty$ and a recursion for all moments of $O_n$ for finite $n$. This variable gives an upper bound for the minimal clade size [2], which is not observable in real data. In applications to genetics, it has been shown to be useful to lower classification errors in genealogical model selection [10].

Equations

$${\lambda }_{b,k}={\int }_0^1{x}^{k-2}(1-x{)}^{b-k}\Lambda (dx)$$

$${\lambda }_b=\sum _{k=2}^b\left(\genfrac{}{}{0ex}{}{b}{k}\right){\lambda }_{b,k}.$$

$$T_n^{(i)}\stackrel{d}{=}{E}_n^{(i)}+M^{(i)}$$

$$L^{(i)}_{n}:=\max\{k,\mbox{ jump $K_{n,i}(k)$ happens earlier than $T^{(i)}_{n}$}\}\in[\kappa_{n,i}],$$

$$\mathcal{C}_{n,i}[L^{(i)}_{n}]=\{j\in[n]:\mbox{$j$ shares all non-private mutations of $i$}\}$$

$$O_n(i):=|{\mathcal{C}}_{n,i}[L_n^{(i)}]|.$$

$$\frac{O_n}{n}\stackrel{a.s.}{\to }S$$

$$E(S^k)=1-\sum _{r=2}^{k+1}{a}_{k+1,r}\frac{\frac{\theta }{2}}{{\lambda }_r+\frac{\theta }{2}},$$

$$E(S)=\frac{\Lambda ([0,1])}{\Lambda ([0,1])+\frac{\theta }{2}},E(S^2)=1-\frac{3}{2}\frac{\frac{\theta }{2}}{\Lambda ([0,1])+\frac{\theta }{2}}+\frac{1}{2}\frac{\frac{\theta }{2}}{{\lambda }_3+\frac{\theta }{2}}.$$

$$\underset{n\to \infty }{\lim }\frac{O_n}{n}=\underset{n\to \infty }{\lim }\frac{{\mathcal{C}}_{n,1}(M)}{n}=f_1(M)(=f_1(M^{(1)}))\text{mbox}a.s.,$$

$$E(S^k)={\int }_0^{\infty }E((f_1(t){)}^k)\frac{\theta }{2}{e}^{-\frac{\theta }{2}t}dt.$$

$$E(S^k)=1-\sum _{r=2}^{k+1}{a}_{k+1,r}{\int }_0^{\infty }{e}^{-{\lambda }_{r}t}\frac{\theta }{2}{e}^{-\frac{\theta }{2}t}dt=1-\sum _{r=2}^{k+1}{a}_{k+1,r}\frac{\frac{\theta }{2}}{{\lambda }_r+\frac{\theta }{2}},$$

$$\frac{O_n}{n}\stackrel{a.s.}{\to }S,$$

$$S\stackrel{a.s.}{=}\sum _{k\in \mathbb{N}}{B}_k{J}_k.$$

$$P\stackrel{d}{=}Gamma({\theta }^{\prime },1),$$

$$S\stackrel{d}{=}\frac{{\sum }_{p\in {\mathcal{P}}_1}p}{{\sum }_{p\in {\mathcal{P}}_1}p+{\sum }_{p\in {\mathcal{P}}_0}p}=:\frac{X}{X+Y}.$$

$$P(g(M_1^{\prime })\le t)=P(M_1^{\prime }\le g^{-1}(t))=1-e^{-\frac{\theta }{2}{g}^{-1}(t)},$$

$$E(S^k)=$$

$$=$$

$$E(S^k)=1-\frac{\theta }{2}\sum _{r=2}^{k+1}{a}_{k+1,r}{\int }_0^{\infty }(1+\rho t{)}^{-\frac{\theta }{2\rho }-1}{e}^{-\left(\genfrac{}{}{0ex}{}{r}{2}\right)t}.$$

$$\frac{O_n}{n}\stackrel{a.s.}{\to }S$$

$$E(S)$$

$$=1-\sum _{k\in \mathbb{N}}{\left(1-\frac{\Lambda ([0,1])}{{\mu }_{-1}}\right)}^k{\left(\frac{{\mu }_{-1}}{{\mu }_{-1}+\frac{\theta }{2}}\right)}^{k-1}\frac{\frac{\theta }{2}}{{\mu }_{-1}+\frac{\theta }{2}}$$

$$=1-\frac{\theta }{2{\mu }_{-1}}\sum _{k\in \mathbb{N}}{\left(\underset{=:a}{\underbrace{\left(1-\frac{\Lambda ([0,1])}{{\mu }_{-1}}\right)\left(\frac{\frac{\theta }{2}}{\frac{\theta }{2}+{\mu }_{-1}}\right)}}\right)}^k$$

$$=1-\frac{\theta }{2{\mu }_{-1}}\frac{a}{1-a}$$

$$E(X_n^j)=\frac{\frac{\theta }{2}}{\frac{\theta }{2}+{\lambda }_n}+\frac{1}{\frac{\theta }{2}+{\lambda }_n}\sum _{k=2}^n\left(\genfrac{}{}{0ex}{}{n}{k}\right){\lambda }_{n,k}E(X_{n-k+1}^j+\sum _{i=1}^j{a}_{i-1,j,k}\frac{X_{n-k+1}^i}{n-k+1}).$$

$$E(O_n^j)=\sum _{k=2}^n\left(\genfrac{}{}{0ex}{}{n-1}{k-1}\right)\frac{{\lambda }_{n,k}}{{\lambda }_n}\sum _{i=0}^j(a_{i,j.k}E(X_{n-k+1}^i)+\frac{(n-k+1)1_{i=j}+b_{i,j,k}}{n-k}E(O_{n-k+1}^i)),$$

$$E(X_n^j)$$

$$=$$

$$E(O_n^j)$$

$$=\sum _{k=2}^n\frac{\left(\genfrac{}{}{0ex}{}{n-1}{k-1}\right){\lambda }_{n,k}}{{\lambda }_n}E((k-1+X_{n-k+1}{)}^j)$$

$$+\sum _{k=2}^n\frac{\left(\genfrac{}{}{0ex}{}{n-1}{k}\right){\lambda }_{n,k}}{{\lambda }_n}E(O_{n-k+1}^j(1-1_{C_n}))$$

$$+\sum _{k=2}^n\frac{\left(\genfrac{}{}{0ex}{}{n-1}{k}\right){\lambda }_{n,k}}{{\lambda }_n}E((k-1+O_{n-k+1}{)}^j{1}_{C_n})$$

Full text

1. Introduction

The potential for adaption of organisms to diverse environments is based on their genetic diversity. Moreover, the specific historic pattern of adaptation and demography leaves distinct marks in the genetic diversity of a sample taken from a population of said organisms. When observing the genetic diversity of a single non-recombining part of the genome, the diversity can be described by the inheritance pattern of the mutations on the genealogical tree of the sample, usually given by a Poisson point process on the genealogy. Modelling the genealogy is thus an important aspect of modelling genetic diversity. Usually, the exact genealogy is not known and cannot be reconstructed perfectly from observed genetic data (an example is provided later). Thus, genealogy models are usually defined as random variables on the set of possible genealogical trees.

Here, we are concerned with the genealogical tree of $n$ alleles, i.e. the genetic information of a sample of size $n$ from a genomic region. Coalescent theory provides a rich class of genealogical tree models for a sample of $n$ alleles as much as elegant tools and a convenient setting up for statistical inference. In particular, Kingman’s coalescent [13] and the larger family of coalescents with multiple collisions [18, 20] were widely studied in the past decades. This class of Markov processes on the set of partitions of $[n]:=\{1,\ldots,n\}$ is characterized by a finite measure $\Lambda$ on $[0,1]$, justifying their name of $\Lambda$-$n$-coalescents. If a $\Lambda$-$n$-coalescent has $b$ blocks, any given $k$-tuple of them will merge at rate

[TABLE]

and the rate for the next coalescence event is

[TABLE]

The starting partition is $\{1\},\ldots,\{n\}$. The genealogical tree is recovered from the partition-valued process by first starting $n$ branches from the leaves $1,\ldots,n$. Then, any merger of partition blocks corresponds to a joining of branches in a node (ancestor), where a single new branch starts. Partition blocks correspond to branches in the genealogical tree, the time a partition block is not merged gives the length of this branch. We refer to [11] for a survey.

Genealogical trees of the alleles in a genetic region come with an interpretation of relatedness and genetic similarity: An allele $i$ is more closely related to allele $j$ than to allele $l$ if the common ancestor of $(i,j)$ appears more recently than the common ancestor of $(i,l)$, while the path lengths between leaves (alleles) measure the time available to accumulate mutations that decrease genetic similarity. Several statistics aim to capture these aspects and their biological meaning. For instance, the minimal clade size of an allele $i$ gives the number of closest relatives of $i$, see [2]. Another example is the length of the external branch of $i$, i.e. the waiting time for the first merger of $i$, which gives a measure of the genetic uniqueness of an individual [19]. The mathematical properties of the minimal clade size [2, 9, 23, 8], the length of an external branch [2, 3, 5], as much as the family (partition block) sizes at this time [2, 24], have been analysed recently. In these works, asymptotic and exact behaviours are obtained for various examples of exchangeable coalescents.

However, these statistics cannot be observed directly from the genetic data. We will illustrate this for the minimal clade. By a clade we denote the set of all alleles that share a specific ancestor, and the minimal clade of $i$ is the smallest clade including $i$. Assume the infinite-sites model of mutation, each mutation causes a change at a different position in the genomic region. Further assume that we know the ancestral state at the genomic region, i.e. we can identify mutations as changes compared to the ancestral state. Any clade can only be observed if there is at least one mutation that all its members share. This mutation is inherited from the common ancestor, thus has to be placed on the branch that connects this ancestor further towards the root of the genealogy (the most recent common ancestor of the whole sample). Thus, we can only observe a clade if there is a mutation on the branch directly above of the ancestor defining this clade.

Instead of looking at the minimal clade of an allele $i$, one could consider the smallest clade which includes $i$ that can be observed from the data. We considered the sizes of these clades for all alleles sampled, the minimal observable clade sizes, in [10]. There we could show that they provide an additional set of statistics that faciliates the inference of a well-fitting genealogy model when coupled with standard statistics of genetic diversity as the site frequency spectrum.

In this article, we study some mathematical properties of the minimal observable clade size for an individual $i$. Its asymptotic behaviour for any $\Lambda$-$n$-coalescent for $n\to\infty$ as well as a recursion for all moments for finite $n$ are established. For the Bolthausen-Sznitman coalescent, which provides a somewhat universal genealogical model for populations under strong selection, see e.g. [17], [4], [22], we can show that the minimal observable clade size is asymptotically Beta-distributed.

2. A formal definition of the minimal observable clade size

Let $[n]:=\{1,\ldots,n\}$ for $n\in\mathbb{N}$, $[n]_{0}=[n]\cup\{0\}$. For any $\Lambda$-$n$-coalescent and a sampled allele $i\in[n]$, define

• •

$\mathcal{C}_{n,i}(t)$ as the partition block $i$ is in at time $t$ (a size-biased pick of a block of the $n$-coalescent at time $t$)

• •

$\kappa_{n}$ as the total number of jumps of the $\Lambda$-$n$-coalescent, $\kappa_{n,i}$ as the total number of jumps (the block of) $i$ participates in

• •

$K_{n,i}(0)(=0),K_{n,i}(1),\ldots,K_{n,i}(\kappa_{n,i})(=\kappa_{n})\in[\kappa_{n}]_{0}$ as the successive indices of jumps in the $\Lambda$-$n$-coalescent in which the block of $i$ is involved

• •

$\mathcal{C}_{n,i}[k]$ as the partition block $i$ is in at the time of its $k$th jump $K_{n,i}(k)$, $k\in[\kappa_{n,i}]_{0}$. $\mathcal{C}_{n,i}[0]=\{i\}$, $\mathcal{C}_{n,i}[1]$ is the minimal clade of $i$, $\mathcal{C}_{n,i}[\kappa_{n,i}]=[n]$.

Given the $\Lambda$-$n$-coalescent tree, we set mutations on its branches via a homogeneous Poisson point process with rate $\frac{\theta}{2}$. Mutations are interpreted under the infinite sites model, each mutation hits a site not hit by any other mutation, producing a new type. The new type is called derived type in contrast to the ancestral type of the most recent ancestor of the sample. Mutations on external branches are affecting only one individual, we will call these private mutations; they can also be referred to as singleton mutations. All other mutations are called non-private mutations. Since we are interested in the mutations carried by individual $i$, we have to record the mutations from $t=0$ to the time back to the most recent common ancestor of the sample (the root of the genealogy) on the path of $i$. Let $T^{(i)}_{n}$ be the waiting time until the first (youngest) mutation on the path of $i$ that is non-private, i.e. does not fall on the external branch which ends in $i$ (which has length $E^{(i)}_{n}$). If we continue the path of $i$ after reaching the most recent common ancestor as a single ancestral line indefinitely (which we will do from now on), we have

[TABLE]

for an independent exponential random variable $M^{(i)}$ with rate $\frac{\theta}{2}$. Let

[TABLE]

i.e. the $L^{(i)}_{n}$th jump that $i$ participates in is the last jump of it before $T^{(i)}_{n}$. The minimal observable clade of $i$ is then given by

[TABLE]

The definitions are equivalent since all non-private mutations of $i$ are inherited from the youngest ancestor of $i$ that bears at least one mutation on the branch connecting it to the next older ancestor. If $i$ has no non-private mutations, $\mathcal{C}_{n,i}[L^{(i)}_{n}]=[n]$ almost surely, since in this case $T^{(i)}_{n}$ is larger than the time back to the most recent common ancestor.

The statistic we are interested in is the size of the minimal observable clade of an allele $i$

[TABLE]

See Figure 1 for an example. Due to exchangeability, the distribution of $O_{n}(i)$ does not depend on $i$, we can even choose $i$ randomly without changing the distribution. For ease of notation, we fix the allele we are interested in to allele 1 and abbreviate $O_{n}:=O_{n}(1)$.

Remark 2.1*.*

Since the partition block including $i$ can only increase in size over time, the minimal clade $\mathcal{C}_{n,i}[1]$ is a subset of the minimal observable clade $\mathcal{C}_{n,i}[L^{(i)}_{n}]$ for $i\in[n]$. Thus, the size $M_{n}(i)$ of the minimal clade of $i\in[n]$ satisfies $M_{n}(i)\leq O_{n}(i)$. See Figure 1 for an example.

3. Asymptotics of the observable clade size

Asymptotically for sample size $n\to\infty$, the probabilistic structure of $O_{n}$ simplifies considerably. First, we focus on coalescents without dust, a class which includes Beta($2-\alpha$,$\alpha$)-$n$-coalescents for $\alpha\in(1,2)$ [21], Kingman’s $n$-coalescent ($\Lambda=\delta_{0}$) and the Bolthausen-Sznitman $n$-coalescent ($\Lambda$ uniform on $[0,1]$). A $\Lambda$-$n$-coalescent has dust if and only if $\mu_{-1}:=\int_{[0,1]}x^{-1}\Lambda(dx)<\infty$, see [18].

Theorem 3.1**.**

Let $O_{n}$ be defined for $\Lambda$-$n$-coalescents such that $\mu_{-1}=\infty$ (without dust) and with mutation rate $\frac{\theta}{2}$. We have

[TABLE]

for $n\to\infty$ with $S>0$ a.s.. $S$ is distributed as the size of the block of individual 1 (alternatively a size-biased pick of a block size) at a random time $M\stackrel{{\scriptstyle d}}{{=}}Exp\left(\frac{\theta}{2}\right)$. The distribution of $S$ is uniquely determined by its moments

[TABLE]

where $\lambda_{r}$ is the total rate of the $\Lambda$-coalescent in a state with $r$ blocks and $a_{k+1,r}$ is a rational function of $\lambda_{2},\ldots,\lambda_{k+1}$, defined as in [18, Prop. 29]. In particular,

[TABLE]

Proof.

Let $E^{(1)}_{n}$ be the waiting time for the first collision of individual 1. By the consistency of the $n$-coalescents, we have $E^{(1)}_{2}\geq E^{(1)}_{3}\geq\ldots$. By a slight adaptation of [18, Prop. 26], we see that $E^{(1)}_{n}\stackrel{{\scriptstyle d}}{{\to}}0$ for $n\to\infty$. Since $(E^{(1)}_{n})_{n\in\mathbb{N}}$ is monotonically decreasing, this convergence also holds almost surely.

All mutations of individual 1 on any $n$-coalescent lie on the path of leaf 1 to the root of the coalescent tree and are consistent for different values of $n$ (any mutation on the path to the root in the $m$-coalescent is also a mutation on this path in every $n$-coalescent with $n>m$), since for $m<n$, the $m$-coalescent (seen as a tree) is the subtree of the $n$-coalescent which is spanned by the leaves $[m]$, including mutations. Thus, we can represent the mutations of individual 1 on the $n$-coalescent by one common homogeneous Poisson process for all $n$, independent of the $n-$coalescents, on $[0,\infty)$ with rate $\frac{\theta}{2}$ This gives a new representation for $T_{n}^{(1)}$, it is the smallest Poisson point $T$ with $T\geq E^{(1)}_{n}$. Let $M$ be the smallest Poisson point overall. Since $E^{(1)}_{n}\to 0$ a.s. for $n\to\infty$, for any realisation of the coalescent there exists a $n_{0}\in\mathbb{N}$ s.t. $E^{(1)}_{m}<M$ a.s. for all $n\geq n_{0}$. This shows that $T_{n}^{(1)}=M$ a.s. for $m\geq n_{0}$, which implies

[TABLE]

where $f_{1}(t)$ is the (asymptotic) frequency of the block individual 1 is in at time $t\geq 0$ in the $\Lambda$-coalescent (with values in $\mathbb{N}$, see [18]). The existence of the limit follows from Kingman’s correspondence, since the coalescent stopped at the random time $M$ (independent of the coalescent) gives an exchangeable partition of $\mathbb{N}$. Since we have a coalescent without dust, we have no singleton blocks a.s. at any time $t>0$ and a (potentially infinite) number of blocks with a.s. positive frequencies, again due to Kingman’s correspondence. This shows $f_{1}(M)>0$ a.s.. Due to exchangeability, the distribution of $f_{1}(M)$ is the same as if we would make a size-biased pick from all blocks present.

Denoting $f_{1}(M)$ by $S$, consider the moments

[TABLE]

From [18, Eq. (50) and Prop. 29] we see that $E((f_{1}(t))^{k})=1-\sum_{r=2}^{k+1}a_{k+1,r}e^{-\lambda_{r}t}$ due to a connection to the exchangeable partition function of the coalescent at time $t$. Thus, we have

[TABLE]

which is Eq. (4). Using the explicit values $a_{2,2}=1,a_{3,2}=\frac{3}{2},a_{3,3}=-\frac{1}{2}$ (essentially from [18, Eq. (39),(40)]) and $\lambda_{2}=\Lambda([0,1])$ yields the first two moments. Since $S$ takes values in $[0,1]$, its distribution is uniquely determined by its moments. ∎

For the special case of the Bolthausen-Sznitman coalescent, the law of $S$ can be identified

Theorem 3.2**.**

For the Bolthausen-Sznitman $n$-coalescent,

[TABLE]

for $n\to\infty$ where $S\stackrel{{\scriptstyle d}}{{=}}Beta\left(\frac{1}{1+\frac{\theta}{2}},\frac{\frac{\theta}{2}}{1+\frac{\theta}{2}}\right)$.

Proof.

[18, Corrolary 16] shows that for the Bolthausen-Sznitman coalescent, $(f_{1}(t))_{t\geq 0}$ jumps at independent standard exponential times with ranked jump sizes given by a Poisson-Dirichlet distribution with parameters $(0,1)$. The set of jump times is independent of the set of jump sizes. Comparing this with Eq. (5), we see that to compute $S\stackrel{{\scriptstyle a.s.}}{{=}}f_{1}(M^{({1})})$, we need to sum the sizes of all jumps of $f_{1}$ that happen before or at $M^{({1})}$. Consider the jump times $(T_{k})_{k\in\mathbb{N}}$ of $f_{1}$ ordered according to the rank of their jump sizes $(J_{k})_{k\in\mathbb{N}}$. Define $B_{k}:=1_{\{T_{k}\leq M^{({1})}\}}$. Hence, we have $P(B_{k}=1)=(1+\frac{\theta}{2})^{-1}$. We can now express

[TABLE]

In other words, $S$ can be seen as summing up a random thinning of a standard Poisson-Dirichlet distributed random variable.

The random variable $S$ is Beta distributed. To see this, we will use the construction of the $PD(0,\theta^{\prime})$ distribution from [14], which is also summarised in [1, Section 4.11]. Consider the points $\mathcal{P}:=(P_{k})_{k\in\mathbb{N}}$ of a Poisson point process on $[0,\infty)$ with mean measure $\nu_{\theta^{\prime}}:=\frac{\theta^{\prime}e^{-x}}{x}dx$. Then, the size-ordered and normalized points $\left(\frac{P_{[k]}}{P}\right)_{k\in\mathbb{N}}$ with $P=\sum_{i\in\mathbb{N}}P_{k}$ have the Poisson-Dirichlet distribution and are independent of

[TABLE]

where $Gamma(\alpha,\beta)$ is the Gamma distribution with shape parameter $\alpha$ and rate $\beta$.

We choose $\theta^{\prime}=1$ and make the correspondence between the ranked and normalised points $(P_{k})_{k\in\mathbb{N}}$ and the jump sizes $(J_{k})_{k\in\mathbb{N}}$. To express Eq. (7), we give each point $P_{k}$ a mark $B_{k}\in\{0,1\}$. Marks are independent from $(P_{k})_{k\in\mathbb{N}}$ and from one another. We set the probability to be marked to $m:=P(B_{k}=1)=(1+\frac{\theta}{2})^{-1}$ for all $k\in\mathbb{N}$. $(P_{k},B_{k})_{k\in\mathbb{N}}$ is a marked Poisson process. The Colouring Theorem [15, Section 5.1] shows that all points $P_{k}$ with marks 1 form a Poisson point process $\mathcal{P}_{1}$ with mean measure $m\nu_{1}$, while all points with mark 0 form a Poisson point process $\mathcal{P}_{0}$ with $(1-m)\nu_{1}$.

We can now alternatively express (7) as

[TABLE]

where $X$ and $Y$ are independent due to the independence of $\mathcal{P}_{0}$ and $\mathcal{P}_{1}$. Since the mean measures of $\mathcal{P}_{0}$ and $\mathcal{P}_{1}$ are of the form $\frac{\theta^{\prime}e^{-x}}{x}dx$ with $\theta^{\prime}$ equal to $m$ and $1-m$, Eq. (8) yields $X\stackrel{{\scriptstyle d}}{{=}}Gamma(m,1)$ and $Y\stackrel{{\scriptstyle d}}{{=}}Gamma(1-m,1)$. Thus, $\frac{X}{X+Y}$ is Beta-distributed with parameters $m$ and $1-m$. ∎

Remark 3.3*.*

Theorem 3.1 can be generalised for some time-changed $\Lambda$-$n$-coalescents without dust, which appear when modelling genealogies in Cannings models with moderate fluctuations in population size, see [12], [16], [25] and [7]. Let $(\Pi^{(n)}_{g(t)})_{t\geq 0}$ be a time-changed $\Lambda$-$n$-coalescent, where $g(t):=\int_{0}^{t}\mu(s)ds$ with continuous $\mu:[0,+\infty)\mapsto[0,+\infty)$, which includes some time changes proposed for $\Lambda$-$n$-coalescents in the references above. Observe that $g$ is continuous, monotone and invertible with differentiable inverse. The time-changed $\Lambda$-$n$-coalescent is still exchangeable. The almost sure convergence of $n^{-1}O_{n}$ for the time-changed $\Lambda$-$n$-coalescent works analogously as in Theorem 3.1. The time of the first merger of individual 1 is $g(E^{(1)}_{n})$, thus also converges to 0 almost surely. The waiting time $M^{\prime}_{1}$ for the first mutation on the path of 1 to the root is an $Exp(\frac{\theta}{2})$-distributed random variable, but on the time-changed path of $(\Pi^{(n)}_{g(t)})_{t\geq 0}$. Thus, the limit of $n^{-1}O_{n}$ for the time-changed $\Lambda$-$n$-coalescent is the frequency of the block containing 1 at time $M^{\prime}_{1}$ in $(\Pi^{(n)}_{g(t)})_{t\geq 0}$. This can also be expressed as $f_{1}(g(M^{\prime}_{1}))$, where $f_{1}$ is said frequency in the $\Lambda$-$n$-coalescent $(\Pi^{(n)}_{t})_{t\geq 0}$. The distribution of $g(M^{\prime}_{1})$ is given by

[TABLE]

which has density $t\mapsto\frac{\theta}{2}\mu(g^{-1}(t))^{-1}e^{-\frac{\theta}{2}g^{-1}(t)}$. Analogously to (6) we can thus express, in terms of the $a_{k,r}$ from Theorem 3.1, the $k$th moment of $n^{-1}O_{n}$ for the $n$-coalescent with exponential growth as

[TABLE]

As an example, consider Kingman’s $n$-coalescent with exponential growth with rate $\rho=0$. From [12], we see that $\mu(t)=e^{\rho t}$ and thus $g^{-1}(t)=\rho^{-1}\log(1+\rho t)$. This leads to moments

[TABLE]

Now, consider coalescents $(\Pi^{(n)}_{t})_{t\geq 0}$ with dust which stay infinite. An example for this are Dirac $n$-coalescents with $\Lambda=\delta_{p}$, $p\in(0,1)$ [6]. For $\Lambda$-coalescents with dust, staying infinite is equivalent to $\Lambda(\{1\})=0$.

Theorem 3.4**.**

Let $O_{n}$ be defined for $\Lambda$-$n$-coalescents with $\mu_{-1}<\infty$ (with dust), $P(\limsup_{n\to\infty}|\Pi^{(n)}_{t}|=\infty\ \forall\ t>0)=1$ and with mutation rate $\frac{\theta}{2}$. We have

[TABLE]

for $n\to\infty$ with $S>0$ a.s.. We have $E(S)=1-\frac{\theta}{2\mu_{-1}}\frac{a}{1-a}$ with $a=\left(1-\frac{\Lambda([0,1])}{\mu_{-1}}\right)\left(\frac{\frac{\theta}{2}}{\frac{\theta}{2}+\mu_{-1}}\right)$.

Proof.

From [8, Thm. 1], we see that the asymptotic frequency of the block of 1 forms an increasing jump-hold process $f_{1}$ with $f_{1}\stackrel{{\scriptstyle a.s.}}{{=}}\lim_{n\to\infty}\frac{\mathcal{C}_{n,1}(t)}{n}$ with values in $[0,1]$, positive jumps and with i.i.d. $Exp(\mu_{-1})$ waiting times between jumps. It fulfills $E(f_{1}[k])=1-(1-\frac{\Lambda([0,1])}{\mu_{-1}})^{k}$, where $f_{1}[k]$ is the value of $f_{1}$ at its $k$th jump. We record between which indices of jumps $K$, $K+1\in\mathbb{N}$ of $f_{1}$ the waiting time $T^{(1)}_{n}$ for the first non-private mutation falls. From the proof of [8, Cor. 1], we know that there exists $n_{0}\in\mathbb{N}$ so that $E^{(1)}_{n}$ equals the time of the first jump of $f_{1}$ for all $n\geq n_{0}$ almost surely. Similarly to the proof of Theorem 3.1, we just have to trace back the first mutation after this first jump whose time of appearance does not depend on $n$. This implies that $T^{(1)}_{n}$ falls between the same $K,K+1$ for all $n\geq n_{0}$. Thus we have $\lim_{n\to\infty}n^{-1}O_{n}=f_{1}[K]$ a.s., where $f_{1}[K]$ is the state of $f_{1}$ at the $K$th jump. We only need to find the distribution of $K$. For $n\geq n_{0}$, $T^{(1)}_{n}$ is the waiting time for the first jump of $f_{1}$ plus an independent $Exp(\frac{\theta}{2})$ random variable $M^{(1)}$. Using that the waiting times $(T_{1,k})_{k\in\mathbb{N}}$ between the jumps of $f_{1}$ are i.i.d., we have $K=1+Y$, where $Y$ is defined by $\sum_{i=1}^{Y}T_{1,i+1}\leq M<\sum_{i=1}^{Y+1}T_{1,i+1}$ and thus $Y\stackrel{{\scriptstyle d}}{{=}}Geo\left(\frac{\frac{\theta}{2}}{\mu_{-1}+\frac{\theta}{2}}\right)$ on $\mathbb{N}_{0}$. This yields $K\stackrel{{\scriptstyle d}}{{=}}Geo\left(\frac{\frac{\theta}{2}}{\mu_{-1}+\frac{\theta}{2}}\right)$ on $\mathbb{N}$. We compute

[TABLE]

∎

4. Recursions for the moments

To obtain recursive formulae for the moments of $O_{n}$ for a $\Lambda$-$n$-coalescent, we first need to introduce $X_{n}$, the size of the block of 1 at the exponential clock $M$ of rate $\frac{\theta}{2}$ in the $n$-coalescent.

Theorem 4.1**.**

Let $j\geq 1$. The $j$th moments of $X_{n}$ and $O_{n}$ satisfy the following recursions: $E(X_{1}^{j})=1,E(O_{2}^{j})=2^{j}$ and

[TABLE]

and

[TABLE]

where $a_{-1,j,k}=0$, $a_{i,j,k}=\binom{j}{i}(k-1)^{j-i},$ and ${b_{i,j,k}=a_{i-1,j,k}-a_{i,j,k}.}$

Proof.

Our proofs rely on tracking the number of blocks involved in the first jump of the $\Lambda$-$n$-coalescent, with some additional condition(s). Let us first prove (9). Let $T_{n}$ be the waiting time for the first coalescence in the $n$-coalescent which is $Exp\left(\lambda_{n}\right)$-distributed. Let $A_{n}$ be the event that the first block merged is part of the block of 1 stopped at $M$.

[TABLE]

Now let us turn to the proof of (10). Let $B_{n}=\{K_{n,1}(1)>1\}$ be the event that 1 does participate in the first coalescence event. Also let $C_{n}$ be the event that the first block merged is part of the observed clade.

[TABLE]

Expanding, we obtain the result. ∎

Remark 4.2*.*

For Kingman’s $n$-coalescent, (9) and (10) considerably simplify. In particular the two first moments of $X_{n}$ are

[TABLE]

and

[TABLE]

and the two first moments of $O_{n}$ are

[TABLE]

and

[TABLE]

FF was funded by DFG grant FR 3633/2-1 through Priority Program 1590: Probabilistic Structures in Evolution.