format: text

text:
# Introduction

## Introduction

The basic mechanisms which create genetic variation in biological evolution are mutation and recombination. They are counteracted by selection, which removes variation. Genetic information may be quite generally described in terms of a collection of linearly ordered sites (i.e. a sequence of sites), each of which is occupied by an element of a given (finite or infinite) set which we denote as site space; if this set is finite, it is often termed alphabet. A specific sequence is also called type.

Mutation is treated as a random state change of a site variable, which occurs independently at every site. Recombination occurs on the occasion of sexual reproduction, and refers to the creation of ‘offspring’ sequences from two (randomly chosen) ‘parental’ ones, where a subset of the ‘maternal’ sites is combined with the complementary set of the ‘paternal’ sites, and the linear ordering along the sequence is maintained. This process is realized through one, or a number of, crossover events, where the two parental strands are interlaced between a pair of neighbouring sites. An important feature of recombination is that it removes dependencies between sites, known as linkage disequilibria in genetics. Finally, selection is caused by the flourishing of fit individuals at the expense of less fit ones.

We consider an infinite population of sequences which evolves under the joint action of mutation, selection or recombination, or of any combination thereof. This is to be considered as the infinite population limit (IPL) of the stochastic process alluded to, and defines a deterministic dynamical system for probability measures (in discrete or in continuous time). It describes the time evolution of the measure with probabilistic certainty, see \[21, Ch. 11\], and Thm. 2.1 of it in particular. Although there are many interesting and important questions connected with finite populations, we focus on the differential equation of the deterministic limit here, which we will call IPL equation from now on. In particular, we will not employ the traditional discrete dynamical systems, but follow the continuous route along the lines of Kimura  and Akin , which happens to be much less developed than it ought to be, see also  for a recent review.

Mutation is a linear process and straightforward to deal with. Selection involves some nonlinearity, which is due to norm conservation under the dynamics, but this nonlinearity may be removed through a simple transformation. Recombination contains a very different source of nonlinearity, which is due to the fact that pairs of objects are involved in the process, and is much harder to treat. Nevertheless, if both state space and time are discrete, a procedure (known as Haldane linearization, see  and \[36, Ch. 6\]) is available which transforms the dynamical system (exactly) into a linear one. It involves a multilinear transformation of the probabilities to a new set of variables, namely certain linkage disequilibria, which describe the deviation from statistical independence of sites. These variables decay independently and geometrically, i.e. they decouple and diagonalize the dynamics. Unfortunately, however, the procedure is cumbersome since it relies on recursions, and no closed form is available for the transformation in the general case.

In a previous paper , the special case of single crossover events was considered, where offspring sequences are composed from one maternal and one paternal segment. This scenario is particularly relevant in molecular evolution, where crossover events are rare, and it is most consistently described in continuous time. For discrete site spaces, and with the help of the corresponding vector space structure, the linearizing transform could be given in closed form with the help of elementary methods from multilinear algebra.

The aim of this article is to further develop this approach in a systematic measure-theoretic setting which also incorporates more general site spaces and does not require an explicit coordinatization. We will essentially start from the deterministic IPL equation and construct its solution explicitly, first for recombination only. The so-called Möbius inversion principle will then give a simple approach to the calculation of a suitable (and, in particular, complete) set of linkage disequilibria. It will then turn out that mutation and even selection may be included in the framework, provided fitness is additive, meaning that the fitness of any type may be decomposed into a sum of independent contributions of its individual sites, i.e. if there is no interaction between sites. Such results may be helpful for the solution of the corresponding inverse problem, i.e. the determination of recombination rates from experimental data, e.g. observed patterns of linkage disequilibria along sequences .

The exposition will be more explicit than needed for a purely mathematical audience, and we also try to give rather precise references to background material we use. We hope that the article will become more self-contained this way and that it is also accessible for readers with a more biological background.

The structure of the paper is as follows. After some preliminaries in Section 1, we will briefly summarize the description of mutation through an IPL equation on the space of positive measures in Section 2, followed by some general remarks on measure-valued IPL equations. The core of the article is Section 3, where we solve, step by step, the IPL equation for recombination and construct an explicit solution of the abstract Cauchy problem, together with a closed form of the corresponding linkage disequilibria. The latter is based on an application of the inclusion-exclusion principle via Möbius inversion (a supplement is given in the Appendix). Section 4 combines mutation and recombination. Section 5 deals with selection and recombination, with some emphasis on the role of mean fitness as a Lyapunov function. Finally, Section 6 ties together all three evolutionary forces — still giving an explicit solution, expressions for the linkage disequilibria, and asymptotic properties. We close with some afterthoughts mainly aimed at the relationship to models in discrete time.

## 1 Preliminaries

If $`X`$ is a locally compact space (by which we always mean to include the Hausdorff property), we use $`_+(X)`$ to denote the set of finite positive regular Borel measures on $`X`$, with $`0_+(X)`$. Likewise, $`(X)`$ is the vector space of real (or signed) finite regular Borel measures. It is a Banach space under the norm $`\omega =|\omega |(X)`$ where $`|\omega |`$ denotes total variation measure. Due to the Riesz-Markov representation theorem, $`(X)`$ can also be viewed as the dual of $`C_{\mathrm{}}(X,)`$, the Banach space of real-valued continuous functions which vanish at infinity, equipped with the usual supremum norm, see \[39, Thm. IV.18\], as well as \[41, Ch. 6\] and \[39, Ch. IV.4\] for general background material. Note that $`(X)`$ with the variation norm $`.`$ is actually a Banach lattice, and this gives access to the highly developed theory of positive operators . We will mainly be interested in the closed convex subsets $`_+^m(X):=\{\omega _+(X)\omega (X)=m\}`$, and in $`𝒫(X)=_+^1(X)`$ in particular, the set of probability measures on $`X`$. Note that, for positive measures $`\omega `$, we simply have $`\omega =\omega (X)`$.

If the Borel $`\sigma `$-algebra of $`X`$ is generated by a family of sets that is closed under finite intersections, a regular Borel measure on $`X`$ is already uniquely specified by its values on the elements of this generating family . This is a property that we will need several times, in particular if $`X=X_1\times X_2`$ is a product space, equipped with the product topology.

###### Fact 1

Let $`\nu ,\nu ^{}`$ be two regular Borel measures on the locally compact product space $`X=X_1\times X_2`$ which coincide on all “rectangles” $`E_1\times E_2`$ where $`E_1`$ and $`E_2`$ each run through the Borel sets of $`X_1`$ and $`X_2`$. Then $`\nu =\nu ^{}`$, i.e. $`\nu (E)=\nu ^{}(E)`$ for all Borel sets $`E`$ of $`X`$.

Proof: In view of the above remark, the only obstacle to cope with is the (non-vacuous!) situation when the $`\sigma `$-algebra generated by the rectangles $`E_1\times E_2`$ is not the full Borel $`\sigma `$-algebra of $`X`$. However, the $`\sigma `$-algebra generated by the rectangles contains all Baire sets $`F`$ of $`X`$, because the Baire sets of $`X`$ possess the required Cartesian product property \[9, Lemma 56.2\], and the Borel sets of $`X_i`$ contain the Baire sets of $`X_i`$. The equality of $`\nu `$ and $`\nu ^{}`$ now follows from \[9, Thm. 62.1\] (this rests on the fact that every Baire measure has a unique extension to a regular Borel measure). $`\mathrm{}`$

Standard examples of locally compact spaces include the compact ones, such as any finite set or the closed interval $`[0,1]`$, but also $`^k`$ and $`^{\mathrm{}}`$ with $`k,\mathrm{}0`$, and arbitrary combinations thereof. These certainly cover all meaningful parameter spaces to be expected in biological applications.

If $`X`$ is a finite set (which is an important case in population genetics), $`𝒫(X)`$ is a simplex. If the cardinality of $`X`$ is $`M`$, this simplex has dimension $`M1`$, i.e. any probability measure is a unique convex linear combination of the $`M`$ extremal measures that constitute the vertices of the simplex. If $`X=\{1,\mathrm{},M\}`$, they are denoted by $`e_i`$, $`i=1,\mathrm{},M`$, and fixed by their values on singleton sets, $`e_i(\{j\})=\delta _{i,j}`$. In other words, any $`\omega 𝒫(X)`$ is of the form $`\omega =_{i=1}^Ma_ie_i`$ with all $`a_i0`$ and $`a_1+\mathrm{}+a_M=1`$. This provides the canonical coordinatization of this situation.

The set (or state space) $`X`$ that we need will have a product structure, described on the basis of sites. For later convenience, we use $`N=\{0,1,\mathrm{},n\}`$ for the set of sites, i.e. we start counting with $`0`$ here. To site $`i`$, we attach the locally compact space $`X_i`$, and our state space is then

$$X=X_0\times X_1\times \mathrm{}\times X_n,$$
(1)
which is still locally compact. One Banach space of measures to show up is the space $`(X)`$ with the corresponding variation norm $`.`$. Note that $`(X)`$ contains the (algebraic) tensor product space $`^{}:=_{i=0}^n(X_i)`$, and also its completion (here, the closure in the given $`.`$-norm of $`(X)`$). To simplify notation, the latter will also be denoted by $`^{}`$, because we shall only deal with Banach spaces here. Recall that $`^{}`$ contains the product measures $`\omega =\omega _0\mathrm{}\omega _n`$ with $`\omega _i(X_i)`$, but also all (finite) linear combinations of measures of this kind. Because we consider the completion, also all measures are contained which can be approximated with such linear combinations in the norm. All probability measures of product form are in this space, but note that the single measures in the product need not be probability measures themselves.

If $`X_i=\{1,\mathrm{},M_i\}`$ is finite, for all $`0in`$, $`X`$ is still a finite set, with $`M=_{i=0}^nM_i`$ elements. Then, $`(X)=^{}`$, and this is simply a real vector space of dimension $`M`$. $`(X)=^{}`$ is also true for $`X`$ discrete. In this case, the action of operators in tensor product form is well defined. In general, if $`^{}(X)`$, one can still go beyond $`^{}`$ under certain circumstances, e.g. by including integrals (rather than finite sums) of product measures. However, we do not want to enter this rather technical discussion, and refer to \[35, Ch. IX.6\] and \[42, Ch. IV.7\] for some background material, and to \[18, Ch. 13\] for some of the problems that are related to these difficulties.

$`X`$ finite is the case most frequently studied in the theory of sequence evolution, and it was the motivation for this work, see  and references given there. However, many results hold in greater generality, which we want to cover in view of potential applications in quantitative genetics. There, the space $`X_i`$ often is a state space such as $``$, or a compact subset thereof. In this case, $`^{}`$ is a true subspace of $`(X)`$, which has to be taken care of later on (occasional restrictions of $`X`$ to a finite set will be mentioned explicitly).

The main reason for using the above set $`N`$ of sites is that we will need ordered partitions of $`N`$, which are uniquely specified by a set of cuts or crossovers. The possible cut positions are at the links between sites, which we denote by half-integers, i.e. by elements of the set $`L=\{\frac{1}{2},\frac{3}{2},\mathrm{},\frac{2n1}{2}\}`$. We will use Latin indices for sites and Greek indices for links, and the implicit rule will always be that $`\alpha =\frac{2i+1}{2}`$ is the link between site $`i`$ and $`i+1`$.

With this notation, the ordered partitions of $`N`$ are in one-to-one correspondence with the subsets of $`L`$ as follows. If $`A=\{\alpha _1,\mathrm{},\alpha _p\}L`$, let $`N_A`$ denote the ordered partition

$$\{0,\mathrm{},\alpha _1\},\{\alpha _1,\mathrm{},\alpha _2\},\mathrm{},\{\alpha _p,\mathrm{},n\}$$

where $`\alpha `$ ($`\alpha `$) denotes the largest integer below $`\alpha `$ (the smallest above $`\alpha `$). In particular, we have $`N_{\mathrm{}}=N`$ and $`N_L=\{\{0\},\mathrm{},\{n\}\}`$. With this definition, it is clear that $`N_B`$ is a refinement of $`N_A`$ if and only if $`AB`$. Consequently, the lattice of ordered partitions of $`N`$ corresponds to the Boolean algebra of the finite set $`L`$, denoted by $`(L)`$, cf. \[1, Ch. I.2\]. We prefer this notation to that with partitions, as it is easier to deal with. If $`AB`$, we will write $`BA`$ for $`BA`$, and $`\stackrel{\mathrm{\_}\mathrm{\_}}{A}`$ for the set $`LA`$.

This setup allows us to use the powerful tool of Möbius inversion from combinatorial theory \[1, Ch. IV.2\], which is a systematic way to employ the inclusion-exclusion principle. If $`f`$ and $`g`$ are mappings from $`(L)`$ to $``$ which are, for all $`AL`$, related by

$$g(A)=\underset{BA}{}f(B),$$
(2)
then this can be solved for $`f`$ via the inversion formula \[1, Thm. 4.18\]

$$f(A)=\underset{BA}{}g(B)\mu (B,A)$$
(3)
with the Möbius function $`\mu (B,A)=(1)^{|AB|}`$, where $`|AB|`$ stands for the cardinality of the set $`AB`$. For $`B`$ not a subset of $`A`$, we set $`\mu (B,A)=0`$ which makes the Möbius function into an element of the so-called incidence algebra, see \[1, Ch. IV.1\] for details. It is important to note that Möbius inversion is not restricted to functions, it also applies to bounded operators.

## 2 Mutation and Markov generator

The description of mutation is rather straight-forward. Let us start from a finite population. Since we are working in continuous time, we assume an independent Poisson clock for each individual member of a (finite) population, and a mutation occurs for an individual whenever its clock rings, according to prescribed mutation rates between (finitely many) types or states. Since the individuals are independent, this is a simple Markov process for each of them. If we now go to the infinite population limit, the time evolution of the probability measure for the types is, almost surely, described by a (deterministic) ordinary differential equation (ODE). This is the so-called IPL equation, compare \[21, Thm. 11.2.1\] for a general justification, which we will also rely on below. For the simple mutation case, this ODE is linear. It clearly coincides with the ODE for the probability measure of the individual Markov process, usually obtained from multiple realizations through the law of large numbers.

Let us consider the case that $`X`$ is a finite state space of cardinality $`|X|=M`$ in more detail, where $`dim_{}((X))=M`$. The mutation rate from state $`\mathrm{}`$ to state $`k`$ is given by $`Q_{k,\mathrm{}}=Q_k\mathrm{}`$, where we already consider $`Q`$ as a mapping acting on the corresponding probabilities, resp. measures. The rate matrix $`Q`$ is a Markov generator, i.e. it has non-negative entries everywhere except on its diagonal, and vanishing column111In contrast to the standard probability literature, we adopt the transposed version here since we are considering the situation from the (linear) operator point of view. sums. The time evolution is then fully described by the Markov semigroup $`\{\mathrm{exp}(tQ)t0\}`$, see \[21, Ch. 1.1 and Ch. 4.2\]. We shall usually assume that $`Q`$ is irreducible, i.e. it is possible to reach every state from any other one. In this case, the equilibrium state is unique and given by the properly normalized 0-eigenvector of the generator $`Q`$. It can actually be given in closed form, see \[24, Lemma 6.3.1\].

If $`X`$ has the product structure introduced above, our mutation process is supposed to be of a more special form, for biological reasons. We assume that mutation happens at all sites in parallel and independently from one another, so that our generator has the form

$$Q=\underset{i=0}{\overset{n}{}}Q_i$$
(4)
where each $`Q_i`$ is, in a properly coordinatized way, the tensor product of a rate matrix at site $`i`$ and unit matrices of matching dimension everywhere else, i.e.

$$Q_i=\mathbf{\hspace{0.33em}1}_{M_0}\mathrm{}\mathrm{𝟏}_{M_{i1}}q_i\mathrm{𝟏}_{M_{i+1}}\mathrm{}\mathrm{𝟏}_{M_n}$$
(5)
where $`q_i`$ is a local rate matrix (of dimension $`M_i`$) for the state space $`X_i`$, acting on $`(X_i)`$. The rate matrices $`Q_i`$ clearly commute with one another. Note also that $`Q`$ of (4) is irreducible if and only if all the $`q_i`$ are. The Markov semigroup inherits the tensor product structure, i.e. we have

$$\mathrm{exp}(tQ)=\underset{i=0}{\overset{n}{}}\mathrm{exp}(tQ_i)=\underset{i=0}{\overset{n}{}}\mathrm{exp}(tq_i).$$
(6)
In view of our following description of recombination, we prefer to avoid an explicit coordinatization here, so we will not use matrix notation. This simply means that we have to reinterpret the generator $`Q`$ as a linear operator on $`(X)`$. Nothing of the above actually changes, we only have to read $`Q`$ (or $`q_i`$) as a linear mapping on $`(X)`$ (or on $`(X_i)`$). The two conditions for $`Q`$ to be a Markov generator now read as follows (the analogous conditions apply to $`q_i`$ in relation to $`(X_i)`$).

1. If $`\nu `$ is a positive measure and $`E`$ any Borel set such that $`\nu (E)=0`$, then $`(Q\nu )(E)0`$.
2. If $`\nu `$ is a positive measure, then $`(Q\nu )(X)=0`$.

The first condition ensures that the semigroup generated by $`Q`$ maps $`_+(X)`$ into itself. Under the present circumstances, where $`Q`$ is bounded and $`(X)`$ is a reflexive Banach space, this condition is necessary and sufficient for the positivity of $`\mathrm{exp}(tQ)`$, see \[5, Thm. 1.11\]. It is sometimes also called the positive minimum principle. The second condition means that the semigroup is Markov, i.e. it preserves the norm of positive measures, and, in particular, maps $`𝒫(X)`$ into itself. In this setting, irreducibility implies that the kernel of the Markov generator $`Q`$ is one-dimensional.

The IPL equation for our simple mutation process222In this linear case, the IPL equation is closely related to the master equation commonly used in the physics literature, see \[29, Ch. 5\] for details. now reads

$$\dot{\omega }=\mathrm{\Phi }_{\mathrm{mut}}(\omega ):=\left(\underset{i=0}{\overset{n}{}}Q_i\right)\omega $$
(7)
which we will take, in generalization of the discrete situation, as the starting point for the analysis of mutation, without tracing it back to an explicit stochastic process. We then obtain, by employing standard results  from the theory of ordinary linear differential equations in (finite-dimensional) Banach spaces (see also Theorem 1 below):

###### Proposition 1

The abstract Cauchy problem of the IPL equation $`(\text{7})`$ with initial condition $`\omega _0𝒫(X)`$ has the unique solution

$$\omega _t=\mathrm{exp}\left(t\underset{i=0}{\overset{n}{}}Q_i\right)\omega _0$$

which is, for $`t0`$, a one-parameter family of probability measures. $`\mathrm{}`$

To formulate a generalization of Prop. 1, let us forget about the product structure for a moment and consider the linear ODE

$$\dot{\omega }=Q\omega $$

with $`Q`$ the generator of a uniformly (or norm) continuous Markov semigroup on $`(X)`$, compare \[19, Ch. I.3\]. This is the case if and only if the linear operator $`Q`$, in addition to satisfying assumptions 1. and 2. from above, is bounded, and hence defined on all of $`(X)`$, see \[19, Cor. II.1.5\]. In particular, we can then write the semigroup in exponential form \[19, Thm. I.3.7\], i.e. as $`\mathrm{exp}(tQ)`$, and the solution as $`\omega _t=\mathrm{exp}(tQ)\omega _0`$. In what follows, we will (non-constructively) assume that a process is given that leads to a bounded generator $`Q`$ which is a linear operator on $`(X)`$, i.e. maps regular Borel measures to regular Borel measures. As long as this is the case, it is sufficient to work with assumptions 1. and 2., even if the space of measures considered is no longer reflexive. The analogue of Prop. 1 then holds on the Banach subspace $`^{}`$, to which we shall restrict our attention whenever $`Q`$ is of the form specified in Eqs. (4) and (5). This makes no difference at all as long as $`X`$ is discrete.

Many results can still be generalized to densely defined generators of strongly continuous semigroups, see \[19, Ch. I.5\], but already the well-posedness of the Cauchy problem needs some thought, compare \[19, Ch. II.6\] for a discussion. Also, the characterization of generators for positive semigroups becomes more involved, see \[5, Ch. 3\]. Usually, one would then rather describe the entire process by means of semigroups on function spaces, compare \[21, Ch. 1.4\]. Since all explicit mutation schemes we have in mind lead to uniformly continuous semigroups, we will not expand on the more general situation.

Let us instead add a few remarks on the general type of IPL equation that arises when recombination and selection are also included. This will also better explain our formulation of mutation, from the point of view of measure-valued differential equations. In what follows, it is sufficient to investigate the first order ODE

$$\dot{\omega }=\mathrm{\Phi }(\omega )$$
(8)
on the Banach space $`(X)`$, where $`\mathrm{\Phi }`$ is a mapping from $`(X)`$ into itself (alternatively, we can study (8) on any closed subspace of $`(X)`$ that is invariant under $`\mathrm{\Phi }`$). Unlike $`\mathrm{\Phi }_{\mathrm{mut}}`$ from (7), $`\mathrm{\Phi }`$ need not be linear, and it is the nonlinear cases below that we are most interested in. The three properties we will meet below are:

* The mapping $`\mathrm{\Phi }`$ is (globally) Lipschitz.
* If $`\nu _+(X)`$, i.e. $`\nu `$ is a positive measure, and $`E`$ any Borel set such that $`\nu (E)=0`$, then we have $`\left(\mathrm{\Phi }(\nu )\right)(E)0`$.
* For any $`\nu _+(X)`$, we have $`\left(\mathrm{\Phi }(\nu )\right)(X)=0`$.

It is clear that our formulation of mutation constitutes a linear example of such a mapping.

###### Theorem 1

If $`\mathrm{\Phi }:(X)(X)`$ satisfies (A1), the abstract Cauchy problem of the ODE $`(\text{8})`$, with initial condition $`\omega _0(X)`$, has a unique solution. If $`\mathrm{\Phi }`$ also satisfies (A2), the cone $`_+(X)`$ of positive measures is invariant under the semiflow for $`t0`$ $`(`$in other words, $`_+(X)`$ is positive invariant$`)`$. Finally, if $`\mathrm{\Phi }`$ also satisfies (A3), the norm of positive measures is preserved in forward time. In particular, the convex set $`𝒫(X)`$ of probability measures is then positive invariant.

Proof: If $`\mathrm{\Phi }`$ is Lipschitz, we can invoke the Picard-Lindelöf Theorem for ODEs on Banach spaces, see \[4, Thm. 7.6\], so existence and uniqueness of the solution of the abstract Cauchy problem are clear.

If $`\mathrm{\Phi }`$ also satisfies (A2), positive invariance of $`_+(X)`$ follows from a continuity argument, see p. 235 and Thm. 16.5 together with Remark 16.6 of  for a proof. If $`\mathrm{\Phi }`$ is linear, (A2) is the so-called positive minimum principle, and our assertion also follows from \[5, Thm. 1.11\], which uses a functional analytic proof.

Finally, assume $`\mathrm{\Phi }`$ satisfies (A1) – (A3). Let $`\omega _0_+^m(X)`$ be the initial condition and denote the corresponding unique solution of (8) by $`\omega _t`$. Then, $`\omega _t_+(X)`$ for all $`t0`$ by the previous argument, so $`\omega _t=\omega _t(X)`$. This implies $`\frac{\mathrm{d}}{\mathrm{d}t}\omega _t=\left(\mathrm{\Phi }(\omega _t)\right)(X)=0`$ by assumption (A3), so $`\omega _t\omega _0=m`$. This proves the assertion. $`\mathrm{}`$

## 3 Recombination

This section deals with the nonlinear IPL equation for recombination, and is the core of our article. We develop the results step by step here. The combination with mutation will then be rather painless, and an addition of selection will be discussed after that.

### 3.1 Recombination on measures

Let $`X,Y`$ be two locally compact spaces with attached measure spaces $`(X)`$ and $`(Y)`$. If $`f:XY`$ is a continuous function and $`\omega (X)`$, then $`f.\omega :=\omega f^1`$ is an element of $`(Y)`$, where $`f^1(y):=\{xXf(x)=y\}`$ means the preimage of $`yY`$ in $`X`$, with obvious extension to $`f^1(B)`$, the preimage of a subset $`BY`$ in $`X`$. Due to the continuity of $`f`$, $`f^1(B)`$ is a Borel set in $`X`$ if $`B`$ is a Borel set in $`Y`$.

Let $`X=X_0\times \mathrm{}\times X_n`$ be as in Section 1, and let, from now on, $`N`$ and $`L`$ always denote the set of sites and links as introduced there. In this section, we can entirely work with the Banach space $`(X)`$, equipped with the variation norm $`.`$. Let $`\pi _i:XX_i`$ be the canonical projection which is continuous. It induces a mapping from $`(X)`$ to $`(X_i)`$ by $`\omega \pi _i.\omega `$, where $`(\pi _i.\omega )(E)=\omega (\pi _i^1(E))`$, for any Borel set $`EX_i`$. By (slight) abuse of notation, we will use the symbol $`\pi _i`$ also for this induced mapping. It is clear that $`\pi _i`$ is linear and maps positive measures to positive measures of the same norm. As such, it is bounded and hence also continuous. In particular, it maps $`𝒫(X)`$ to $`𝒫(X_i)`$ and may then be understood as marginalization. Likewise, we can start from any (ordered) index set $`IN`$ and define a projector $`\pi _I:(X)(X_I)`$ with $`X_I:=\times _{iI}X_i`$. With this notation, $`X_N=X`$. We will frequently also use the abbreviation $`\pi _{<\alpha }`$ for the projector $`\pi _{\{1,\mathrm{},\alpha \}}`$, and $`\pi _{>\alpha }`$ for $`\pi _{\{\alpha ,\mathrm{},n\}}`$. These objects may be understood as ‘cut and forget’ operators, since they give the distribution of what is left after a cut is made at $`\alpha `$, and the trailing resp. leading segment is discarded.

This now enables us to introduce the elementary recombination operator, or recombinator as we will call it from now on, $`R_\alpha :(X)(X)`$, for $`\alpha L`$. If $`\omega =0`$, $`R_\alpha (\omega ):=0`$, and otherwise

$$R_\alpha (\omega ):=\frac{1}{\omega }((\pi _{<\alpha }.\omega )(\pi _{>\alpha }.\omega ))$$
(9)
which is a (partial) product measure. Here and in what follows, we tacitly identify (if necessary) a product measure with its unique extension to a regular Borel measure on $`X`$, which is justified by Fact 1. The following property is now an immediate consequence of the definition.

###### Fact 2

The recombinator $`R_\alpha `$ maps $`_+(X)`$ into itself and preserves the norm of positive measures. In particular, it maps $`𝒫(X)`$ into itself. $`\mathrm{}`$

Let us comment on the choice of (9). Being composed of the cut-and-forget operators for the leading and the trailing ends, $`R_\alpha (\omega )`$ has the interpetation of a ‘cut-and-relink operator’, which describes a cut at $`\alpha `$, followed by (random) reunion of the resulting segments.

At first sight, it might appear more natural to drop the prefactor $`1/\omega `$. However, the norm of a positive measure $`\omega `$ would then not be preserved unless $`\omega =1`$. In view of later extensions, it is more desirable not to be restricted to probability measures, and that is why we prefer (9) which makes $`R_\alpha `$ positive homogeneous of degree 1,

$$R_\alpha (a\omega )=|a|R_\alpha (\omega ),$$
(10)
for arbitrary $`a`$. Note, however, that $`R_\alpha `$ is not a linear operator, not even when restricted to $`_+(X)`$.

###### Fact 3

Let $`\alpha L`$. The recombinator $`R_\alpha `$ satisfies $`R_\alpha (\omega )\omega `$, for all $`\omega (X)`$, and is $`(`$globally$`)`$ Lipschitz on $`(X)`$.

Proof: Let us first observe that, for arbitrary $`\omega ,\omega ^{}(X)`$ and $`\alpha L`$, we obtain the inequality

$$(\pi _{<\alpha }.\omega )(\pi _{>\alpha }.\omega ^{})\omega \omega ^{},$$

which is a simple consequence of Hahn’s decomposition for real measures, see \[41, Thm. 6.14\], applied separately to the factors of the product measure. For $`0\omega (X)`$, we then have

$$R_\alpha (\omega )=\frac{(\pi _{<\alpha }.\omega )(\pi _{>\alpha }.\omega )}{\omega }\omega ,$$

with equality for positive measures, as stated in Fact 2. Clearly, we also have $`R_\alpha (0)=0`$, so that the first assertion follows.

Let $`\omega ,\omega ^{}(X)`$. If one of them is the 0-measure, say $`\omega ^{}=0`$, we have $`R_\alpha (\omega )R_\alpha (\omega ^{})=R_\alpha (\omega )\omega =\omega \omega ^{}`$. So we may assume both $`\omega `$ and $`\omega ^{}`$ to be different from 0 and hence to have positive norm. With the above inequalities, we can now employ the following 3$`\epsilon `$-type argument

$`R_\alpha (\omega )R_\alpha (\omega ^{})`$
$`=`$ $`{\displaystyle \frac{(\pi _{<\alpha }.\omega )(\pi _{>\alpha }.(\omega \omega ^{}))}{\omega }}+{\displaystyle \frac{(\pi _{<\alpha }.(\omega \omega ^{}))(\pi _{>\alpha }.\omega ^{})}{\omega ^{}}}`$
$`+({\displaystyle \frac{1}{\omega }}{\displaystyle \frac{1}{\omega ^{}}})[(\pi _{<\alpha }.\omega )(\pi _{>\alpha }.\omega ^{})]`$
$``$ $`{\displaystyle \frac{(\pi _{<\alpha }.\omega )(\pi _{>\alpha }.(\omega \omega ^{}))}{\omega }}+{\displaystyle \frac{(\pi _{<\alpha }.(\omega \omega ^{}))(\pi _{>\alpha }.\omega ^{})}{\omega ^{}}}`$
$`+|{\displaystyle \frac{1}{\omega }}{\displaystyle \frac{1}{\omega ^{}}}|(\pi _{<\alpha }.\omega )(\pi _{>\alpha }.\omega ^{})`$
$``$ $`2\omega \omega ^{}+\left|\omega \omega ^{}\right|\mathrm{\hspace{0.33em}\hspace{0.33em}3}\omega \omega ^{}.`$
Together, this gives the second assertion, with Lipschitz constant $`3`$. $`\mathrm{}`$

In view of Fact 2, it makes sense to investigate the properties of the recombinators restricted to the positive cone $`_+(X)`$. The crucial property which underlies our later analysis is the following.

###### Proposition 2

The elementary recombinators, when restricted to $`_+(X)`$, are idempotents and commute with one another. In other words, we then have $`R_\alpha ^2=R_\alpha `$ and $`R_\alpha R_\beta =R_\beta R_\alpha `$ for arbitrary $`\alpha ,\beta L`$.

Proof: The statement is trivial for the action on $`\nu =0`$. So, let $`\nu >0`$ be a (strictly) positive measure. We then have $`\nu (X)=\nu `$ and obtain

$`\pi _{<\alpha }.((\pi _{<\alpha }.\nu )(\pi _{>\alpha }.\nu ))`$ $`=`$ $`\nu (\pi _{<\alpha }.\nu )`$
$`\pi _{>\alpha }.((\pi _{<\alpha }.\nu )(\pi _{>\alpha }.\nu ))`$ $`=`$ $`\nu (\pi _{>\alpha }.\nu )`$
in $`_+(X_{<\alpha })`$ resp. $`_+(X_{>\alpha })`$ where we adopt the same index convention for sets as we did for projectors. Using $`R_\alpha (\nu )=\nu `$ from Fact 2 and the linearity of the mappings $`\nu (\pi .\nu )`$, one can now apply the definition of the elementary recombinators to check explicitly that

$$R_\alpha \left(R_\alpha (\nu )\right)=R_\alpha (\nu ).$$

For commutativity, we may again assume $`\nu >0`$ and also $`\alpha <\beta `$. Then

$`\pi _{<\alpha }.((\pi _{<\beta }.\nu )(\pi _{>\beta }.\nu ))`$ $`=`$ $`\nu (\pi _{<\alpha }.\nu )`$
$`\pi _{>\alpha }.((\pi _{<\beta }.\nu )(\pi _{>\beta }.\nu ))`$ $`=`$ $`((\pi _{\{\alpha ,\mathrm{},\beta \}}.\nu )(\pi _{>\beta }.\nu )).`$
The first equation can be verified directly, as in the previous case. The second can easily be checked on Borel sets of the product form $`E=E_{\{\alpha ,\mathrm{},\beta \}}\times E_{>\beta }`$, followed by an application of Fact 1. Combining these intermediate results, one obtains

$$R_\alpha \left(R_\beta (\nu )\right)=\frac{1}{\nu ^2}((\pi _{<\alpha }.\nu )(\pi _{\{\alpha ,\mathrm{},\beta \}}.\nu )(\pi _{>\beta }.\nu ))=R_\beta \left(R_\alpha (\nu )\right),$$

which proves our assertion. $`\mathrm{}`$

Remark: In view of positive homogeneity of the recombinators, see Eq. (10), it would have been sufficient to prove our assertions on $`𝒫(X)`$. The above version, however, shows quite clearly where, and how many, normalization factors $`\nu `$ appear in the tensor products. If we restrict ourselves to probability measures below, one should keep this in mind for extending arguments to the full cone, $`_+(X)`$.

A close inspection of the proof of Proposition 2 shows that we have simultaneously proved the following useful property.

###### Lemma 1

Let $`\nu 𝒫(X)`$ and $`\alpha L`$. For all $`\beta L`$ with $`\beta \alpha `$, we have $`\pi _{<\alpha }.\left(R_\beta (\nu )\right)=\pi _{<\alpha }.\nu `$. Similarly, $`\pi _{>\alpha }.\left(R_\beta (\nu )\right)=\pi _{>\alpha }.\nu `$, for all $`\beta \alpha `$. $`\mathrm{}`$

### 3.2 The IPL equation and its solution

Let us start with a brief description of the recombination process for finite $`X`$, and a population of $`m`$ individuals, each of the form $`x=(x_0,x_1,\mathrm{},x_n)`$ with $`x_iX_i`$. Every individual carries a Poisson clock at each link $`\alpha L`$, with parameters $`\varrho _\alpha >0`$, which do not depend on the individual. If the clock at link $`\alpha `$ of the individual $`x`$ rings, a random partner $`y`$ is picked from the population for recombination at that link. The recombined pair is then $`(x_0,\mathrm{},x_\alpha ,y_\alpha ,\mathrm{},y_n)`$ and $`(y_0,\mathrm{},y_\alpha ,x_\alpha ,\mathrm{},x_n)`$.

To describe the entire population, let $`Z_x(t)`$ be the random variable that gives the number of $`x`$-individuals at time $`t`$, and $`Z(t)`$ the combined random vector with components $`Z_x(t)`$. Hence, if $`Z(t)=z`$, and $`xy`$, we can have transitions from $`z`$ to $`zu_xu_y+u_{(x_{<\alpha },y_{>\alpha })}+u_{(y_{<\alpha },x_{>\alpha })}`$, where we use our short hand notation for indices, and $`u_x`$ to denote the unit vector corresponding to $`x`$. Such a transition occurs at rate $`\varrho _\alpha z_xz_y/(mz_x)`$.

Note that this process implies instant mixing of all (geno-)types in the population. This is an idealization which neglects that maternal and paternal genes stay together for the lifetime of an individual. Nevertheless, this is a good and realistic model if recombination events are rare on the time scale of the individual life span. This is certainly true if our sites belong to the DNA sequence of a single gene, or a few adjacent genes. It is then well justified to describe recombination in terms of these first order effects only.

Let us look at the influence of increasing $`m`$, whence we write $`Z^{(m)}(t)`$ to indicate dependence on system size. As $`m\mathrm{}`$, the sequence of random processes $`Z^{(m)}(t)/m`$ converges almost surely to the solution of a differential equation with initial condition $`Z^{(m)}(0)/m`$ (resp. its limit as $`m\mathrm{}`$), see \[21, Thm. 11.2.1\]. The corresponding IPL equation \[6, Eq. 2.5\], reformulated in our measure-theoretic setting, reads

$$\dot{\omega }=\mathrm{\Phi }_{\mathrm{rec}}(\omega ):=\underset{\alpha L}{}\varrho _\alpha \left(R_\alpha \mathrm{𝟏}\right)(\omega ).$$
(11)
In line with our strategy for the mutation processes, we take this nonlinear ODE as the general starting point for the recombination analysis on product spaces $`X`$ built from arbitrary locally compact spaces $`X_i`$. We will assume that $`\varrho _\alpha >0`$, for all $`\alpha L`$, without loss of generality (if $`\varrho _\alpha =0`$, remove the link at $`\alpha `$, absorb the pair $`(\alpha ,\alpha )`$ into a single site, and identify $`X_\alpha \times X_\alpha `$ with the state space at that site, thus reducing the number of sites (and links) by one).

###### Proposition 3

The abstract Cauchy problem of the IPL equation $`(\text{11})`$ has a unique solution. Furthermore, $`_+(X)`$ is positive invariant under the flow, with the norm of positive measures preserved. In particular, $`𝒫(X)`$ is positive invariant.

Proof: Consider $`\dot{\omega }=\mathrm{\Phi }_{\mathrm{rec}}(\omega )`$, which is a special case of (8), so we want to apply Theorem 1. By Fact 3, $`\mathrm{\Phi }_{\mathrm{rec}}`$ is Lipschitz, so assumption (A1) is satisfied.

Let $`\nu _+(X)`$, i.e. $`\nu (E)0`$ for all Borel sets $`EX`$. Let $`E`$ be any Borel subset of $`X`$ such that $`\nu (E)=0`$. Then

$$\mathrm{\Phi }_{\mathrm{rec}}(\nu )(E)=\underset{\alpha L}{}\varrho _\alpha R_\alpha (\nu )(E)\mathrm{\hspace{0.33em}0}$$

because each $`R_\alpha (\nu )`$ is a positive measure and all $`\varrho _\alpha >0`$ by assumption, so (A2) is satisfied.

Finally, with Fact 2, it is easy to check that $`\mathrm{\Phi }_{\mathrm{rec}}(\nu )(X)=0`$ for any positive measure $`\nu `$, which shows that assumption (A3) is satisfied, too. Theorem 1 then establishes our claims. $`\mathrm{}`$

The difficulty in solving $`(\text{11})`$ stems from the nonlinearity of the right-hand side, so $`\mathrm{\Phi }_{\mathrm{rec}}`$ cannot be considered as the generator of an exponential semigroup. It is, however, rather natural to expect that the solution should still have a rather similar structure, as the $`R_\alpha `$ are at least positive homogeneous of degree one and commute with one another. Let us therefore, for any $`GL`$, introduce the composite recombinators

$$R_G:=\underset{\alpha G}{}R_\alpha .$$
(12)
They are well-defined on $`_+(X)`$ due to Proposition 2, while an order of the product has to be specified otherwise. In any case, $`R_G(\omega )\omega `$ for all $`\omega (X)`$. Note that $`R_{\mathrm{}}=1`$ and $`R_{\{\alpha \}}=R_\alpha `$ in this notation. The composite recombinators are again positive homogeneous of degree one. A simple induction argument based on Proposition 2 gives the following result.

###### Corollary 1

On $`_+(X)`$, the composite recombinators satisfy

$$R_GR_H=R_{GH},$$

for arbitrary $`G,HL`$. Furthermore, each $`R_G`$ maps $`_+(X)`$ into itself and preserves the norm of positive measures. $`\mathrm{}`$

Let us pretend for a moment that the idempotents $`R_\alpha `$ were actually linear operators. In such a case, we would get

$$\mathrm{exp}\left(\varrho _\alpha t(R_\alpha \mathrm{𝟏})\right)=\mathrm{exp}(\varrho _\alpha t)\mathrm{𝟏}+\left(1\mathrm{exp}(\varrho _\alpha t)\right)R_\alpha .$$

Taking the product over such terms for all $`\alpha L`$ and expanding it would formally lead to the sum

$$\underset{GL}{}a_G(t)R_G$$

with the coefficient functions

$$a_G(t)=\left(\underset{\alpha \overline{G}}{}\mathrm{exp}(\varrho _\alpha t)\right)\left(\underset{\beta G}{}\left(1\mathrm{exp}(\varrho _\beta t)\right)\right).$$
(13)
It will have a touch of magic below when we prove that this little “derivation” actually gives the correct answer! After we have established our main result in Theorem 2, we will come back to these coefficients and give them a probabilistic interpretation. This will also motivate why they are a very reasonable guess to start with.

As mentioned before, the elementary recombinators are not linear. Nevertheless, they have a related property on convex combinations. If $`\omega =_{i=1}^ka_i\nu _i`$ is a convex linear combination of positive measures $`\nu _i`$ of equal norm, we get

$$R_\alpha (\omega )=\underset{i=1}{\overset{k}{}}a_iR_\alpha (\nu _i)+$$
(14)
where one can show, by a rather straight-forward calculation which we omit here, that the remainder $``$ is given by

$$=\underset{i<j}{}a_ia_j\frac{\nu _i\nu _j}{\omega }R_\alpha (\nu _i\nu _j).$$

This shows that the recombinators are indeed inherently nonlinear, but also that they might act like linear operators on special convex combinations, namely those for which the remainder vanishes. This is precisely what we need to solve our problem.

###### Proposition 4

Let $`\nu `$ be a positive measure, $`\alpha L`$ and $`a_G(t)`$ the coefficient functions of $`(\text{13})`$. Then, for any fixed $`t0`$, we have

$$R_\alpha \left(\underset{GL}{}a_G(t)R_G(\nu )\right)=\underset{GL}{}a_G(t)R_{G\{\alpha \}}(\nu ).$$

Before we prove this result, we formulate a special property of the coefficient functions first. Observe that, for fixed $`t0`$, $`q_\alpha :=\mathrm{exp}(\varrho _\alpha t)`$ is a number between $`0`$ and $`1`$. It can be interpreted as a probability (namely that link $`\alpha `$ has not been hit until time $`t`$). With this, the coefficients read

$$a_G^L=\underset{\alpha LG}{}q_\alpha \underset{\beta G}{}(1q_\beta )\mathrm{\hspace{0.33em}0}$$
(15)
where we have suppressed the (fixed) time, but added the set of links, $`L`$, as an upper index. We can now formulate a crucial factorization property.

###### Lemma 2

Let $`L=L_1\dot{}L_2`$ be a partition of $`L`$, and set $`G_i=GL_i`$ for an arbitrary $`GL`$. Then, the coefficients of $`(\text{15})`$ satisfy $`a_G^L=a_{G_1}^{L_1}a_{G_2}^{L_2}`$. Furthermore, for any $`L^{^{}}L`$, we have

$$\underset{HL^{^{}}}{}a_H^L^{^{}}=\mathrm{\hspace{0.33em}1}.$$

Proof: Since $`L_1L_2=\mathrm{}`$, the first statement is a direct consequence of the product form of $`a_G^L`$ in Eq. (15). The normalization property can be verified from the probabilistic interpretation mentioned above. If $`1q_\alpha `$ (resp. $`q_\alpha `$) is the probability that link $`\alpha `$ has (resp. has not) been hit, $`a_H^L^{^{}}`$ is the probability that, of the links in $`L^{^{}}`$, precisely $`H`$ is spared. Consequently, $`_{HL^{^{}}}a_H^L^{^{}}`$ is the sum over the probabilities of all possible events, hence equal to 1. Alternatively, this identity can be derived from a simple Möbius inversion argument, as we show below in Fact 5. $`\mathrm{}`$

Proof of Proposition 4: Since the recombinators are positive homogeneous of degree one, it suffices to prove the statement for $`\nu `$ a probability measure. Let $`\alpha L`$ be fixed.

Set $`\omega =_{GL}a_G^LR_G(\nu )`$. Since $`\nu 𝒫(X)`$ implies $`\omega 𝒫(X)`$, we obtain

$`R_\alpha (\omega )`$ $`=`$ $`R_\alpha \left({\displaystyle \underset{GL}{}}a_G^LR_G(\nu )\right)`$
$`=`$ $`(\pi _{<\alpha }.\left({\displaystyle \underset{GL}{}}a_G^LR_G(\nu )\right))(\pi _{>\alpha }.\left({\displaystyle \underset{HL}{}}a_H^LR_H(\nu )\right))`$
$`=`$ $`{\displaystyle \underset{G,HL}{}}a_G^La_H^L((\pi _{<\alpha }.R_G(\nu ))(\pi _{>\alpha }.R_H(\nu )))`$
where we have used the linearity of the mappings $`\pi _{<\alpha }`$ and $`\pi _{>\alpha }`$.

Let us define $`L_1=\{\frac{1}{2},\frac{3}{2},\mathrm{},\alpha \}`$ and $`L_2=LL_1`$, so that $`L=L_1\dot{}L_2`$ is a partition of $`L`$. Also, let $`G_i=GL_i`$ and $`H_i=HL_i`$, for $`G,HL`$. Lemma 1 then tells us that

$`(\pi _{<\alpha }.R_G(\nu ))(\pi _{>\alpha }.R_H(\nu ))`$ $`=`$ $`(\pi _{<\alpha }.R_{G_1}(\nu ))(\pi _{>\alpha }.R_{H_2}(\nu ))`$
$`=`$ $`(\pi _{<\alpha }.R_{G_1H_2}(\nu ))(\pi _{>\alpha }.R_{G_1H_2}(\nu ))`$
$`=`$ $`R_\alpha \left(R_{G_1H_2}(\nu )\right).`$
Inserting this into the previous equation and invoking Lemma 2 repeatedly gives

$`R_\alpha (\omega )`$ $`=`$ $`{\displaystyle \underset{G,HL}{}}a_G^La_H^LR_\alpha \left(R_{G_1H_2}(\nu )\right)`$
$`=`$ $`{\displaystyle \underset{G_1L_1}{}}{\displaystyle \underset{G_2L_2}{}}{\displaystyle \underset{H_1L_1}{}}{\displaystyle \underset{H_2L_2}{}}a_{G_1}^{L_1}a_{G_2}^{L_2}a_{H_1}^{L_1}a_{H_2}^{L_2}R_\alpha \left(R_{G_1H_2}(\nu )\right)`$
$`=`$ $`{\displaystyle \underset{G_1L_1}{}}{\displaystyle \underset{H_2L_2}{}}a_{G_1}^{L_1}a_{H_2}^{L_2}R_\alpha \left(R_{G_1H_2}(\nu )\right)`$
$`=`$ $`{\displaystyle \underset{KL}{}}a_K^LR_\alpha \left(R_K(\nu )\right),`$
which proves our assertion. $`\mathrm{}`$

Remark: Proposition 4 admits the following interpretation. Let $`\nu `$ be a positive measure, with $`\nu =m>0`$. Then, the $`2^{|L|}`$ measures $`R_G(\nu )`$ with $`GL`$ form the vertices of a $`.`$-closed simplex in $`_+^m(X)`$. On some of their convex combinations (in particular along solutions, as we will see shortly), the elementary recombinators $`R_\alpha `$ act linearly. It is this simplex, foliated into solution curves, to which the entire time evolution is constrained, with $`\nu `$ as the initial condition.

The positive measure in Proposition 4 was arbitrary. This means that, when restricting the action of the $`R_\alpha `$’s to $`_+(X)`$, we can formulate the rule on the level of operators. Observe that $`R_\alpha \left(R_K(\nu )\right)=R_{K\{\alpha \}}(\nu )=R_K(\nu ^{})`$ where $`\nu ^{}=R_\alpha (\nu )`$. By a simple induction argument, we thus arrive at

###### Corollary 2

Let $`a_G(t)`$ be the coefficient function of $`(\text{13})`$, and let $`t0`$ be fixed. On $`_+(X)`$, the recombinators satisfy the equation

$$R_H\left(\underset{GL}{}a_G(t)R_G\right)=\underset{GL}{}a_G(t)R_{GH}$$

for arbitrary $`HL`$. $`\mathrm{}`$

We now assume that the initial condition, $`\omega _0`$, is a positive measure and make the following ansatz for the solution of (11):

$$\omega _t=\underset{GL}{}a_G(t)R_G(\omega _0)$$
(16)
with the coefficient functions $`a_G(t)`$ of (13). Note that they do not depend on $`\omega _0`$. The initial values are $`a_{\mathrm{}}(0)=1`$ and $`a_G(0)=0`$ for all $`\mathrm{}GL`$. By Corollary 1, each $`R_G(\omega _0)`$ is a positive measure with the same norm as $`\omega _0`$. This implies that, as long as $`a_G(t)0`$ for $`t0`$, the ansatz for $`\omega _t`$ must form a convex linear combination of positive measures of equal norm if it is a solution of (7). This follows from Eq. (15) together with Lemma 2, or from Fact 5 below.

The time derivative of $`\omega _t`$ of (16) is $`\dot{\omega }_t=_{GL}\dot{a}_G(t)R_G(\omega _0)`$. On the other hand, Proposition 4 means that the $`R_\alpha `$ act linearly on the convex combination (16), and we obtain

$`\mathrm{\Phi }_{\mathrm{rec}}(\omega _t)`$ $`=`$ $`{\displaystyle \underset{\alpha L}{}}\varrho _\alpha (R_\alpha \mathrm{𝟏})(\omega _t)`$
$`=`$ $`{\displaystyle \underset{\alpha L}{}}\varrho _\alpha {\displaystyle \underset{GL}{}}a_G(t)\left(R_{G\{\alpha \}}(\omega _0)R_G(\omega _0)\right)`$
$`=`$ $`{\displaystyle \underset{\alpha L}{}}\varrho _\alpha \left[{\displaystyle \underset{\alpha \text{\$̣G\$}L}{}}a_{G\{\alpha \}}(t)R_G(\omega _0){\displaystyle \underset{\alpha \text{\$̣G\$}L}{}}a_G(t)R_G(\omega _0)\right]`$
$`=`$ $`{\displaystyle \underset{GL}{}}\left[{\displaystyle \underset{\alpha G}{}}\varrho _\alpha a_{G\{\alpha \}}(t){\displaystyle \underset{\beta \overline{G}}{}}\varrho _\beta a_G(t)\right]R_G(\omega _0),`$
where we use the notation $̣G$ in the third step to indicate the summation variable. It is now a straight-forward calculation to check that the coefficients $`a_G(t)`$ of (13) indeed satisfy the equations

$$\dot{a}_G(t)=\underset{\alpha G}{}\varrho _\alpha a_{G\{\alpha \}}(t)\underset{\beta \overline{G}}{}\varrho _\beta a_G(t)$$

and that they constitute a convex combination in (16). Consequently, our ansatz solves the IPL equation (11), and, by Proposition 3, this is the unique solution we are after. We have thus established the following main result.

###### Theorem 2

The ansatz $`(\text{16})`$ solves the IPL equation $`(\text{11})`$ with initial condition $`\omega _0_+(X)`$ if and only if the coefficient functions are given by $`(\text{13})`$, i.e. by

$$a_G(t)=\mathrm{exp}\left(\underset{\alpha \overline{G}}{}\varrho _\alpha t\right)\underset{\beta G}{}\left(1\mathrm{exp}(\varrho _\beta t)\right)$$

for all $`GL`$. $`\mathrm{}`$

Remark: To interpret the coefficient $`a_G(t)`$, let us consider a single individual. Since $`\mathrm{exp}(\varrho _\alpha t)`$ is the probability that link $`\alpha `$ has experienced no crossover event until time $`t`$ (recall that we have assumed a Poisson process of rate $`\varrho _\alpha `$ at link $`\alpha `$), $`a_G(t)`$ may be interpreted as the probability that the set of all links that have, up to time $`t`$, experienced at least one crossover event, is precisely $`G`$.

Note that the above result relies on the assumption of single, independent crossover events, which is described by recombinators that commute. In more general models, with multiple, dependent events, the algebraic structure is rather involved. This makes solutions much more cumbersome, or considerably less explicit in nature (for review, see \[36, Ch. 6\]).

Let us come back to the meaning of Eq. (16) in combination with Theorem 2. If $`\phi _t`$ denotes the flow of the IPL equation (11), we obtain, for all $`t0`$, the identity

$$\phi _t=\underset{GL}{}a_G(t)R_G$$
(17)
which is valid on the cone $`_+(X)`$. As usual, $`\phi _0=\mathrm{𝟏}`$ and $`\phi _t\phi _s=\phi _{t+s}`$, for all $`t,s0`$. This implies the identity

$$a_G(t+s)=\underset{\begin{array}{c}H,KL\\ HK=G\end{array}}{}a_H(t)a_K(s),$$

which can be verified by direct computation. More interestingly, we also have

###### Fact 4

On $`_+(X)`$, the forward flow of $`(\text{11})`$ commutes with the recombinators, i.e. $`R_G\phi _t=\phi _tR_G`$, for all $`t0`$ and $`GL`$.

Proof: Let $`\nu _+(X)`$ and fix $`GL`$. Then

$`R_G\left(\phi _t(\nu )\right)`$ $`=`$ $`R_G\left({\displaystyle \underset{HL}{}}a_H(t)R_H(\nu )\right)={\displaystyle \underset{HL}{}}a_H(t)R_{GH}(\nu )`$
$`=`$ $`{\displaystyle \underset{HL}{}}a_H(t)R_H\left(R_G(\nu )\right)=\phi _t\left(R_G(\nu )\right)`$
by an application of Corollary 2. $`\mathrm{}`$

Once the solution is known, the remaining task is to identify linear combinations of the $`R_H(\omega )`$ that decouple from each other and decay exponentially. To this end, we employ combinatorial techniques to regroup the terms of the solution according to their exponential damping factors. Let us first expand the expression for $`a_G(t)`$,

$$a_G(t)=\underset{KG}{}(1)^{|GK|}\mathrm{exp}\left(\underset{\alpha \overline{K}}{}\varrho _\alpha t\right).$$

This suggests to define new functions $`b_K(t)`$ via

$$b_K(t)=\mathrm{exp}\left(\underset{\alpha \overline{K}}{}\varrho _\alpha t\right),$$
(18)
with the usual convention that the empty sum is 0. In particular, we have $`b_{\mathrm{}}(t)=a_{\mathrm{}}(t)=\mathrm{exp}\left(_{\alpha L}\varrho _\alpha t\right)`$ and $`b_K(0)=1`$ for all $`KL`$. Now, the Möbius inversion of (2) and (3), used backwards, gives us the relation

$$b_K(t)=\underset{GK}{}a_G(t).$$

One immediate consequence is

$$\underset{GL}{}a_G(t)=b_L(t)=\mathrm{\hspace{0.33em}1}.$$
(19)
So, together with the observation that the functions $`a_G(t)`$ of Theorem 2 are always non-negative, we have independently confirmed

###### Fact 5

If $`\omega _0_+(X)`$, the coefficient functions $`a_G(t)`$ of Theorem 2 constitute a convex linear combination of positive measures in Eq. $`(\text{16})`$. $`\mathrm{}`$

The significance of the new functions becomes clear by realizing that there is an analogue on the level of operators. To this end, we rewrite the composite recombinators in terms of new operators via $`R_H=_{\text{\$̣G\$}H}T_G`$ and obtain, by an obvious variant of Möbius inversion,

$$T_G:=\underset{\text{\$̣H\$}G}{}(1)^{|HG|}R_H.$$
(20)
A straight-forward calculation then reveals that

$$\omega _t=\underset{GL}{}a_G(t)R_G(\omega _0)=\underset{KL}{}b_K(t)T_K(\omega _0).$$
(21)
Note that, as a consequence of Eqs. (10), (12) and (20), the operators $`T_G`$ are positive homogeneous of degree one, i.e.

$$T_G(\alpha \nu )=|\alpha |T_G(\nu ).$$
(22)
Let us now introduce new measures $`\nu _G(t):=b_G(t)T_G(\omega _0)`$, which are elements of $`(X)`$, but no longer positive in general.

###### Proposition 5

The signed measures $`\nu _G(t)`$ solve the Cauchy problem

$$\dot{\nu }_G(t)=\left(_{\alpha \overline{G}}\text{ }\varrho _\alpha \right)\nu _G(t)$$

with initial condition $`\nu _G(0)=T_G(\omega _0)`$, for all $`GL`$.

Proof: The result is a direct consequence of the fact that the coefficient functions $`b_G(t)`$ solve the ordinary initial value problems

$$\dot{b}_G(t)=\left(_{\alpha \overline{G}}\varrho _\alpha \right)b_G(t)$$

with initial conditions $`b_G(0)=1`$, see above. $`\mathrm{}`$

So, the transformation (20) resulted in regrouping the terms of the solution to the IPL equation (11) according to their exponential decay factors in time. In particular,

$$\nu _L(t)T_L(\omega _0)=R_L(\omega _0)=\underset{i=0}{\overset{n}{}}(\pi _i.\omega _0)$$

is the unique limit measure of the process starting from $`\omega _0`$. Due to the action of $`R_L`$, it is a complete product measure and reflects total independence, and we obtain $`\omega _tR_L(\omega _0)`$ as $`t\mathrm{}`$ in the $`.`$-topology. This is so because

$$\omega _t\nu _L(t)=\underset{KL}{}b_K(t)T_K(\omega _0)\underset{KL}{}b_K(t)T_K(\omega _0)$$

where all remaining coefficient functions $`b_K(t)`$, i.e. those with $`KL`$, decay exponentially (recall that $`\varrho _\alpha >0`$ for all $`\alpha L`$).

### 3.3 Linkage disequilibria

Starting from the measures $`\nu _G(t)`$, we will now identify a minimal, complete set of variables by evaluating certain $`k`$-point cylinder functions (called $`k`$-point functions from now on) or correlation functions known as linkage disequilibria in genetics. They are important for data analysis because they allow to evaluate associations between sites up to a given order from measured type frequencies, and average over all others by marginalization. This way, a certain amount of stochasticity, which is present in all real (finite) populations, is smoothed out.

Various different definitions of linkage disequilibria are available in the literature (see \[11, p. 183–186\] for an overview). But only special choices decouple (see ), and these are the linkage disequilibria we are after. In view of the applications, we will now restrict ourselves to the case that $`X`$ is a finite set, although the results hold, with only minor modifications, also more generally. Eq. (21) and Proposition 5 suggest to employ the signed measures $`T_G(\omega _0)`$. The corresponding functions $`b_G(t)`$ will then describe their evolution in time.

Let $`j_1,\mathrm{},j_k`$, with $`j_1<\mathrm{}<j_k`$, symbolically denote a cylinder set in $`X=X_N`$ which is specified at sites $`j_i`$, for $`1ik`$. More specifically, these are sets of the product form

$$j_1,\mathrm{},j_k=X_{\{0,\mathrm{},j_11\}}\times \{x_{j_1}\}\times \left[\mathrm{}\right]\times \{x_{j_k}\}\times X_{\{j_k+1,\mathrm{},n\}}$$

where $`[\mathrm{}]`$ contains factors $`\{x_i\}`$ or $`X_i`$ depending on whether $`i`$ appears in $`j_1,\mathrm{},j_k`$ or not. For $`\nu 𝒫(X)`$ and arbitrary $`\alpha L`$, we then have

$`R_\alpha (\nu )\left(j_1,\mathrm{},j_k\right)`$ $`=`$ $`((\pi _{<\alpha }.\nu )(\pi _{>\alpha }.\nu ))\left(j_1,\mathrm{},j_k\right)`$
$`=`$ $`\{\begin{array}{cc}\nu \left(j_1,\mathrm{},j_k\right)\hfill & \text{if }\alpha <j_1\text{ or }\alpha >j_k\hfill \\ \nu \left(j_1,\mathrm{},j_s\right)\nu \left(j_{s+1},\mathrm{},j_k\right)\hfill & \text{if }j_s<\alpha <j_{s+1}\text{.}\hfill \end{array}`$
For later convenience, we also define $`\mathrm{}=X`$ so that $`R_\alpha (\nu )\left(\mathrm{}\right)=1`$.

###### Lemma 3

If $`\nu 𝒫(X)`$, we have $`T_G(\nu )\left(j_1,\mathrm{},j_k\right)=0`$ whenever the set $`\overline{G}`$ contains an element that is less than $`j_1`$ or larger than $`j_k`$.

Proof: Let $`I=\{\beta j_1<\beta <j_k\}`$. Assume there is an $`\alpha \overline{GI}`$. Then

$`T_G(\nu )\left(j_1,\mathrm{},j_k\right)`$
$`=`$ $`{\displaystyle \underset{\text{\$̣H\$}G}{}}(1)^{|HG|}R_H(\nu )\left(j_1,\mathrm{},j_k\right)`$
$`=`$ $`{\displaystyle \underset{\text{\$̣H\$}G}{}}(1)^{|HG|}R_{HI}(\nu )\left(j_1,\mathrm{},j_k\right)`$
$`=`$ $`{\displaystyle \underset{\alpha \text{\$̣H\$}G}{}}\left((1)^{|HG|}+(1)^{|(H\{\alpha \})G|}\right)R_{HI}(\nu )\left(j_1,\mathrm{},j_k\right)`$
where the previous calculation was used in the second step, and summation is over $`H`$. Clearly, the last expression vanishes. $`\mathrm{}`$

Let us now define the time-dependent $`k`$-point functions as

$$F_G^{(t)}(j_1,\mathrm{},j_k)=T_G(\omega _t)\left(j_1,\mathrm{},j_k\right)$$
(23)
for arbitrary $`GL`$, where the notation is again symbolic in that we only specify the positions $`j_i`$, but not the corresponding values. To relate this to Eq. (21), we show

###### Proposition 6

If $`\omega _0_+(X)`$, we have $`T_G(\omega _t)=b_G(t)T_G(\omega _0)`$, for all $`GL`$ and $`t0`$.

Proof: Since $`b_G(0)=1`$, equality holds for $`t=0`$, and the claim follows if we show that $`T_G(\omega _t)`$ and $`b_G(t)T_G(\omega _0)`$ satisfy the same differential equation. With $`\omega _t=\phi _t(\omega _0)`$, compare (17), we obtain

$`{\displaystyle \frac{d}{dt}}T_G(\omega _t)`$ $`=`$ $`{\displaystyle \frac{d}{dt}}{\displaystyle \underset{\text{\$̣H\$}G}{}}(1)^{|HG|}R_H\left(\phi _t(\omega _0)\right)`$
$`=`$ $`{\displaystyle \underset{\text{\$̣H\$}G}{}}(1)^{|HG|}{\displaystyle \frac{d}{dt}}\phi _t\left(R_H(\omega _0)\right)\text{(by Fact }\text{4}\text{)}`$
$`=`$ $`{\displaystyle \underset{\text{\$̣H\$}G}{}}(1)^{|HG|}\mathrm{\Phi }_{\mathrm{rec}}\left(\phi _t(R_H(\omega _0))\right)\text{(by Eq. (}\text{11}\text{))}`$
$`=`$ $`{\displaystyle \underset{\text{\$̣H\$}G}{}}(1)^{|HG|}{\displaystyle \underset{\alpha \overline{G}}{}}\varrho _\alpha \left(R_{H\{\alpha \}}R_H\right)(\omega _t)\text{(by Fact }\text{4}\text{)}`$
$`=`$ $`\left({\displaystyle \underset{\alpha \overline{G}}{}}\varrho _\alpha \right)T_G(\omega _t).`$
The last step is correct because

$$\underset{\alpha \overline{G}}{}\varrho _\alpha \underset{\text{\$̣H\$}G}{}(1)^{|HG|}R_{H\{\alpha \}}=\mathrm{\hspace{0.33em}0}$$

by an argument analogous to the one used in the last step of the proof of Lemma 3. Now, a comparison with Proposition 5 establishes the claim. $`\mathrm{}`$

Even after Lemma 3, there are still too many functions around. It is thus reasonable to select an independent set from them. To see how to do this, assume that we have an index $`\alpha GI`$, with $`I=\{\beta j_1<\beta <j_k\}`$ for a cylinder set of type $`j_1,\mathrm{},j_k`$ as above. Let $`H`$ be a subset of $`L`$ that contains $`G`$, so $`\alpha H`$ in particular, and $`\nu 𝒫(X)`$. Then $`R_H(\nu )=R_\alpha (\nu _H)`$ with $`\nu _H=R_{H\{\alpha \}}(\nu )`$. The little calculation before Lemma 3 now tells us that

$`R_H(\nu )\left(j_1,\mathrm{},j_k\right)`$ $`=`$ $`R_\alpha (\nu _H)\left(j_1,\mathrm{},j_k\right)`$
$`=`$ $`\left[\nu _H\left(j_1,\mathrm{},j_s\right)\right]\left[\nu _H\left(j_{s+1},\mathrm{},j_k\right)\right]`$
$`=`$ $`\left[R_\alpha (\nu _H)\left(j_1,\mathrm{},j_s\right)\right]\left[R_\alpha (\nu _H)\left(j_{s+1},\mathrm{},j_k\right)\right]`$
$`=`$ $`\left[R_H(\nu )\left(j_1,\mathrm{},j_s\right)\right]\left[R_H(\nu )\left(j_{s+1},\mathrm{},j_k\right)\right]`$
where $`j_s<\alpha <j_{s+1}`$. Consequently, defining $`I_1=\{\beta j_1<\beta <j_s\}`$ and $`I_2=\{\beta j_{s+1}<\beta <j_k\}`$, and referring back to (20), we also get

$`T_G(\nu )\left(j_1,\mathrm{},j_k\right)`$
$`=`$ $`{\displaystyle \underset{\text{\$̣H\$}G}{}}(1)^{|HG|}\left[R_H(\nu )\left(j_1,\mathrm{},j_s\right)\right]\left[R_H(\nu )\left(j_{s+1},\mathrm{},j_k\right)\right]`$
$`=`$ $`{\displaystyle \underset{\text{\$̣H\$}G}{}}(1)^{|HG|}{\displaystyle \underset{\begin{array}{c}\text{\$̣K!\$}_1H\\ \text{\$̣K!\$}_2H\end{array}}{}}\left[T_{K_1}(\nu )\left(j_1,\mathrm{},j_s\right)\right]\left[T_{K_2}(\nu )\left(j_{s+1},\mathrm{},j_k\right)\right]`$
$`=`$ $`{\displaystyle \underset{\text{\$̣H\$}G}{}}(1)^{|HG|}{\displaystyle \underset{\begin{array}{c}\text{\$̣K!\$}_1(H\overline{I}_1)\\ \text{\$̣K!\$}_2(H\overline{I}_2)\end{array}}{}}\left[T_{K_1}(\nu )\left(j_1,\mathrm{},j_s\right)\right]\left[T_{K_2}(\nu )\left(j_{s+1},\mathrm{},j_k\right)\right]`$
where Lemma 3 was used in the last step to remove terms that vanish. This equation means that $`T_G(\nu )\left(j_1,\mathrm{},j_k\right)`$, whenever an $`\alpha GI`$ exists, either vanishes (if Lemma 3 applies) or is a polynomial expression in $`\mathrm{}`$-point functions with $`\mathrm{}<k`$.

In the above calculation, $`\nu `$ is an arbitrary probability measure, wherefore the equations apply to $`\omega _t`$, for an arbitrary $`t0`$. Whenever $`GI\mathrm{}`$, the time-dependent $`k`$-point functions are polynomially dependent of $`\mathrm{}`$-point functions with $`\mathrm{}<k`$. Consequently, they do not contain new information. So far, we have:

###### Proposition 7

The $`k`$-point function $`F_G^{(t)}(j_1,\mathrm{},j_k)=T_G(\omega _t)\left(j_1,\mathrm{},j_k\right)`$ can only be non-vanishing and $`(`$polynomially$`)`$ independent from $`\mathrm{}`$-point functions with $`\mathrm{}<k`$ if $`G=\overline{I}=\{\beta <j_1\}\{\beta >j_k\}`$. $`\mathrm{}`$

We choose this collection of $`k`$-point functions as our linkage equilibria.

Let us finally observe that the summation of a $`k`$-point function over all possible values $`x_i`$ of one of the specified $`X_i`$ (i.e. marginalization) reduces it to a $`(k1)`$-point function, so we have one extra (linear) relation. This means that only $`M_i1`$ possible values can be prescribed independently at site $`i`$. On the other hand, given $`j_1,\mathrm{},j_k`$, there is only one way to choose $`G`$ due to Proposition 7, and then there are $`(M_{j_1}1)\mathrm{}(M_{j_k}1)`$ different and independent choices to specify the actual values at the sites. Summing up all these possibilities results in

$$\underset{DN}{}\underset{iD}{}(M_i1)=\underset{i=0}{\overset{n}{}}\left(1+(M_i1)\right)=\underset{i=0}{\overset{n}{}}M_i=|X|.$$

This means that we have singled out the right number of functions. In view of Proposition 7, for $`t`$ arbitrary but fixed, they completely determine the value of the signed measures $`T_G(\omega _t)`$ on all cylinder sets. These, in turn, are closed under finite intersections and generate the full $`\sigma `$-algebra of the (finite) space $`X`$, so all measures $`T_G(\omega _t)`$, and hence also $`\omega _t`$, are uniquely specified, and we have achieved our goal. An explicit example has been worked out in Section 4 of , where the $`k`$-point functions $`F_G^{(t)}(j_1,\mathrm{},j_k)`$ appear as the components of the vector $`z`$ of linkage disequilibria, up to a change of basis in the local site spaces.

If $`X`$ is not a finite set, one has to use a generating family of Borel cylinder sets instead of just singleton sets, and invoke Fact 1. Although there is no simple counting argument, the general structure is still similar.

At this point, one could still argue that $`k`$-point functions w.r.t. the selection of sites, as our $`F_G^{(t)}(j_1,\mathrm{},j_k)`$, should be replaced by proper $`k`$-point correlation functions because these separate off all contributions of functions of lower order, i.e. of $`\mathrm{}`$-point functions with $`\mathrm{}<k`$. This is just another application of the Möbius inversion principle, but one where all partitions (rather than only ordered ones) are needed. We provide the corresponding formulas in the Appendix. If one performs the necessary calculations, one quickly realizes that our previous inclusion-exclusion process w.r.t. ordered partitions of the links has far reaching consequences: most of the potential correction terms simply vanish, as a result of Lemma 3. In particular, we obtain

###### Theorem 3

Let $`S=\{j_1,\mathrm{},j_k\}`$ be a set of site indices, in increasing order and without gaps, and let $`G=\{\alpha <j_1\}\{\alpha >j_k\}`$. Then, the $`k`$-point function $`F_G^{(t)}(j_1,\mathrm{},j_k)=T_G(\omega _t)\left(j_1,\mathrm{},j_k\right)`$ coincides with the corresponding $`k`$-point correlation function as given in Eq. $`(\text{43})`$ of the Appendix.

These functions, for all possible choices of the set $`S`$, form a polynomially independent set of linkage disequilibria.

Proof: We apply Lemma 3 with $`\nu =\omega _t`$. Due to the assumption on $`S`$ versus $`G`$, the right-hand side of Eq. (43) boils down to the one term we already have, because all other terms vanish. Propositions 7 and 6 ensure the polynomial independence of these objects, which are our linkage disequilibria. $`\mathrm{}`$

This result does not extend to all $`k`$-point functions. If, for a given $`k`$-point function, a non-vanishing correction term occurs in the corresponding correlation function, this will, in general, not decay with the same exponential rate as the original $`k`$-point function. So, grouping according to decay rates and according to correlation structures simultaneously is not possible in general. It is a rather remarkable fact that the set of linkage disequlibria is a set of exceptions, and one (as we demonstrated above for the case of discrete state spaces) that completely determines the probability measure.

## 4 Mutation and recombination

In this section, we will just combine the results of the previous two sections. This is possible because, as we will see, mutation and recombination are independent in our approach, i.e. the corresponding operators in the IPL equation commute. This is to be expected given the fact that mutation acts on the sites while recombination works via the links. However, to be able to formulate this in a more general situation than $`X`$ finite or discrete, we now restrict ourselves to the Banach space $`^{}=_{iN}(X_i)`$ which, as explained earlier, is meant as the completion of the algebraic tensor product. In general, it is a (true) Banach subspace of $`(X)`$. Our IPL equation now reads

$$\dot{\omega }=\left(\underset{iN}{}\mu _iQ_i+\underset{\alpha L}{}\varrho _\alpha \left(R_\alpha \mathrm{𝟏}\right)\right)(\omega )$$
(24)
where we have taken the liberty to introduce mutation rates $`\mu _i`$, all of which are assumed to be strictly positive. The idea behind this is to use some standardized version for the mutation operators $`Q_i`$ of (5) so that the $`\mu _i`$ serve as relative coefficients, in line with the usual practice in the biological literature. The linear operators $`Q_i`$ are supposed to be bounded, hence continuous, and thus possess a unique extension to $`^{}`$, compare \[45, Thms. II.1.2 and II.1.5\]. To show consistency, we observe

###### Lemma 4

The Banach space $`^{}`$ is invariant under $`R_\alpha `$, for all $`\alpha L`$, and hence positive invariant under the flow of $`(\text{24})`$.

Proof: It is clear that $`R_\alpha `$ maps a finite linear combination of product measures onto another linear combination of this kind, compare the proof of Prop. 4. Since such linear combinations are dense in $`^{}`$ and $`R_\alpha `$ is Lipschitz on $`(X)`$, it maps the closed subspace $`^{}`$ of $`(X)`$ into itself. The statement on positive invariance is a direct consequence of \[4, Thm. 16.5 and Remark 16.6\]. $`\mathrm{}`$

Let $`𝒫^{}`$ be the subspace of probability measures in $`^{}`$. Referring back to Propositions 1 and 3 and to Theorem 1, the following result is immediate.

###### Proposition 8

The abstract Cauchy problem of the IPL equation $`(\text{24})`$, with initial condition $`\omega _0^{}`$, has a unique solution. The cone $`_+^{}`$ is positive invariant, and the norm of a positive measure is preserved in forward time. In particular, the convex set $`𝒫^{}`$ is positive invariant. $`\mathrm{}`$

To continue, let us call a positive linear operator $`W`$ on $`^{}`$ strictly positive if $`\omega _+^{}`$ with $`\omega >0`$ implies $`W\omega >0`$. The key observation is now

###### Lemma 5

Let $`W`$ be a strictly positive bounded linear operator on $`^{}`$ which has a complete tensor product structure, i.e. $`W=w_0\mathrm{}w_n`$. On $`_+^{}`$, the elementary recombinator $`R_\alpha `$ then commutes with $`W`$, i.e. $`WR_\alpha =R_\alpha W`$. In particular, this is true if $`W=\mathrm{exp}(tQ_i)`$ is an element of a Markov semigroup, as in Section 2, for any $`t0`$, $`iN`$ and $`\alpha L`$.

Proof: Let us first consider the case that $`W`$ preserves the norm of a positive measure $`\nu `$, i.e. $`W\nu =\nu `$. This is also true of $`R_\alpha `$, $`\alpha L`$. Since $`W`$ is linear and $`R_\alpha `$ positive homogeneous of degree 1, it is sufficient to prove the claim on $`𝒫^{}`$. So, let $`\nu 𝒫^{}`$. $`W`$ has a complete tensor product structure, so $`W=W_{<\alpha }W_{>\alpha }`$ in particular. Observe first that $`W_{<\alpha }\pi _{<\alpha }=\pi _{<\alpha }W`$ and $`W_{>\alpha }\pi _{>\alpha }=\pi _{>\alpha }W`$. These relations certainly hold when applied to a product measure $`\nu =\nu _{<\alpha }\nu _{>\alpha }`$, but, due to linearity of all mappings involved here, also on arbitrary (finite) linear combinations of measures of this kind. The latter are dense in $`^{}`$, so that continuity of the mappings establishes the relations, compare \[45, Thm. II.1.5\].

As a consequence, we obtain

$$W((\pi _{<\alpha }.\nu )(\pi _{>\alpha }.\nu ))=(\pi _{<\alpha }.(W\nu ))(\pi _{>\alpha }.(W\nu )),$$

which proves the assertion for the case that $`W`$ preserves the norm of $`\nu `$.

Let us now consider the general case. The proof so far only required that $`W`$ preserved the norm of the single $`\nu `$ under consideration. We employ again positive homogeneity of $`R_\alpha `$. If $`\nu >0`$, we have $`W\nu >0`$ by assumption, so that $`a:=\frac{\nu }{W\nu }>0`$ is well defined. So we obtain $`\nu =aW\nu `$ and

$`WR_\alpha (\nu )`$ $`=`$ $`\frac{1}{a}(aW)R_\alpha (\nu )`$
$`=`$ $`\frac{1}{a}R_\alpha (aW\nu )\text{(by above argument)}`$
$`=`$ $`R_\alpha (W\nu )\text{(by Eq. (}\text{10}\text{))}`$
which proves the first assertion.

The second claim is obvious because elements of a Markov semigroup are strictly positive and because the generators $`Q_i`$, compare Eq. (5), have the required product structure. $`\mathrm{}`$

We can now put together our previous efforts. The obvious form of the solution of (24) is now

$$\omega _t=\mathrm{exp}\left(tQ\right)\underset{GL}{}a_G(t)R_G(\omega _0)$$
(25)
with $`Q=_{i=0}^n\mu _iQ_i`$ and the coefficient functions $`a_G(t)`$ of Theorem 2. The verification that this indeed solves the IPL equation is a simple application of the product rule. Let $`\nu _t:=_{GL}a_G(t)R_G(\omega _0)`$, so that $`\omega _t=\mathrm{exp}(tQ)\nu _t`$. Then we have

$`\dot{\omega }_t`$ $`=`$ $`Q\omega _t+\mathrm{exp}\left(tQ\right)\dot{\nu }_t`$
$`=`$ $`Q\omega _t+\mathrm{exp}\left(tQ\right){\displaystyle \underset{\alpha L}{}}\varrho _\alpha \left(R_\alpha \mathrm{𝟏}\right)(\nu _t)\text{(by Theorem }\text{2}\text{)}`$
$`=`$ $`Q\omega _t+\left({\displaystyle \underset{\alpha L}{}}\varrho _\alpha (R_\alpha \mathrm{𝟏})\right)(\mathrm{exp}(tQ)\nu _t)\text{(by Lemma }\text{5}\text{)}`$
$`=`$ $`\left({\displaystyle \underset{iN}{}}\mu _iQ_i+{\displaystyle \underset{\alpha L}{}}\varrho _\alpha \left(R_\alpha \mathrm{𝟏}\right)\right)(\omega _t).`$
So, together with Proposition 8, we have established:

###### Theorem 4

The unique solution of the IPL equation $`(\text{24})`$, with initial condition $`\omega _0_+^{}`$, is given by $`\omega _t`$ of $`(\text{25})`$, with the coefficient functions $`a_G(t)`$ of Theorem 2. $`\mathrm{}`$

Let us take a closer look at the asymptotic behaviour. Since $`a_G(t)`$ decreases, as $`t\mathrm{}`$, exponentially to $`0`$ unless $`G=L`$, we obtain

$`\omega _t`$ $``$ $`\mathrm{exp}(tQ)\left(a_L(t)R_L(\omega _0)\right)\mathrm{exp}(tQ){\displaystyle \underset{i=0}{\overset{n}{}}}(\pi _i.\omega _0)`$
$`=`$ $`{\displaystyle \underset{i=0}{\overset{n}{}}}\left(\mathrm{exp}(tq_i)(\pi _i.\omega _0)\right)`$
where all neglected terms are of lower order in that they vanish exponentially (recall that $`\mathrm{exp}(tQ)`$ and $`\mathrm{exp}(tq_i)`$ are Markov). This shows that the stationary measure, for any initial measure $`\omega _0`$, is again a complete product measure333Convergence to product measures is also known from various interacting particle systems, compare .. Whether or not there is a unique global equilibrium measure then depends on the properties of the local mutation operators $`q_i`$. In the case that $`X`$ is finite, uniqueness follows if all these generators are irreducible.

What remains to be done is to extend the Möbius trick and to evaluate the linkage disequilibria also for this case. Due to Lemma 5, we can equivalently write $`\omega _t`$ of (25) as

$$\omega _t=\underset{GL}{}a_G(t)R_G\left(\mathrm{exp}(tQ)\omega _0\right).$$
(26)
At any fixed instant of time, $`\mathrm{exp}(tQ)\omega _0`$ is a positive measure, and we can employ Eq. (21) to obtain

$$\omega _t=\underset{KL}{}b_K(t)T_K\left(\mathrm{exp}(tQ)\omega _0\right)$$
(27)
with the functions $`b_K(t)`$ introduced in (18).

If we now assume again that $`X`$ is finite, we can use the $`k`$-point cylinder functions as before to select a finite set of linkage disequilibria that completely determine the solution $`\omega _t`$. They are the functions

$$F_G^t(j_1,\mathrm{},j_k)=T_G(\omega _t)\left(j_1,\mathrm{},j_k\right)$$
(28)
for $`GL`$ and selected cylinder sets $`j_1,\mathrm{},j_k`$ exactly as before.

Since mutation and recombination are independent of each other and the time evolutions commute, we can separate the time decay due to the two processes. The effect is as follows. Recombination is sensitive to sites selected in the cylinder sets, but not to the actual values prescribed there. Mutation, in turn, has a tensor product structure with respect to the sites (which expresses the independence of individual events).

If $`\mathrm{exp}(tQ)`$ is Markov (so that Lemma 5 applies), it is easy to derive (in analogy with the proof of Proposition 6) that

$$\frac{d}{dt}T_G(\omega _t)=\left(Q\underset{\alpha \overline{G}}{}\varrho _\alpha \right)T_G(\omega _t).$$
(29)
This shows how the recombination rates and the eigenvalues of $`Q`$ together determine the fine structure of exponential decay. Note that diagonalizing $`Q`$ (if at all possible) now corresponds to taking appropriate linear combinations of $`T_G(\omega _t)\left(j_1,\mathrm{},j_k\right)`$ for fixed $`G`$ and $`j_1,\mathrm{},j_k`$, but different values prescribed at the sites. For finite $`X`$, this has been worked out in , along with explicit examples.

## 5 Selection

Let us first look at selection in a slightly more general way, i.e. via an IPL equation on $`(X)`$ without explicit reference to its tensor product structure. Let $`P:(X)(X)`$ be a bounded linear operator which generates a positive semigroup. According to \[5, Thm. 1.11\], the latter is true if and only if $`P`$ satisfies our assumption (A2), the positive minimum principle. Consider now the ODE

$$\dot{\omega }=\mathrm{\Phi }_{\mathrm{sel}}(\omega ):=P\omega \frac{P\omega (X)}{\omega }\omega $$
(30)
where $`\mathrm{\Phi }_{\mathrm{sel}}(0)=0`$ is the proper extension of $`\mathrm{\Phi }_{\mathrm{sel}}(\omega )`$ to $`\omega =0`$. This is motivated by the standard selection model (cf. ), where, in properly coordinatized form as indicated in Section 1, $`P`$ is a diagonal matrix which keeps track of the ‘fitness’ of the various states, and $`\frac{P\omega (X)}{\omega }`$ is the ‘mean fitness’ of the population. This model also arises in the infinite population limit of the well-known Moran model, see  or \[22, Ch. 3\]. Here, in a population of $`m`$ individuals with finite state space $`X`$ as described in Section 3, every individual of type $`x`$ reproduces at rate $`r_x`$, and the offspring replaces a randomly chosen individual in the population (possibly its own parent). Therefore, a transition from population state $`z`$ to $`z+u_xu_y`$ occurs at rate $`r_xz_xz_y/m`$. Along the lines of Section 3, the limit $`m\mathrm{}`$ yields a special case of the differential equation (30), where $`P`$ is the diagonal matrix with elements $`r_x`$.

The more general form used here does not only cover more general $`X`$, but also interaction between mutation and reproduction (as opposed to the independent processes considered so far), e.g. the production of mutated offspring on the occasion of reproduction. In any case, the subtraction of the second term on the right hand side of (30) comes from the preservation of total mass, or, in more technical terms, is designed so that $`\mathrm{\Phi }_{\mathrm{sel}}`$ satisfies assumptions (A2) and (A3) from Section 2.

So far, our selection equation seems to imply that selection acts on haploids (i.e. individuals with only one copy of the genetic information per cell). If, however, individuals have two copies that are equivalent and do not interact (the diploid case without dominance), Eq. (30) is replaced by

$$\dot{\omega }=\frac{M(\omega P\omega +P\omega \omega )}{\omega }\frac{\left(M(\omega P\omega +P\omega \omega )\right)(X)}{\omega ^2}\omega ,$$
(31)
where $`M(\mu \nu ):=\nu (X)\mu `$ denotes marginalization with respect to the second factor. In this formulation, the mean fitness is

$$\frac{1}{\omega ^2}\left(M(\omega P\omega +P\omega \omega )\right)(X)=\mathrm{\hspace{0.33em}2}\frac{\omega (X)P\omega (X)}{\omega ^2}.$$

For positive $`\omega `$, the right-hand side of (31) becomes

$$\frac{P\omega (X)}{\omega }\omega +P\omega \mathrm{\hspace{0.17em}2}\frac{P\omega (X)}{\omega }\omega =P\omega \frac{P\omega (X)}{\omega }\omega ,$$

that is, the diploid equation reduces to the haploid one in this case, in the sense that the flow is the same on $`_+(X)`$.

Let us now take a closer look at the differential equation (30).

###### Fact 6

The mapping $`\mathrm{\Phi }_{\mathrm{sel}}:(X)(X)`$ is $`(`$globally$`)`$ Lipschitz.

Proof: Consider $`\omega ,\omega ^{}(X)`$. If one of them is the zero measure, $`\omega ^{}`$ say, we get

$`\mathrm{\Phi }_{\mathrm{sel}}(\omega )\mathrm{\Phi }_{\mathrm{sel}}(0)`$ $`=`$ $`\mathrm{\Phi }_{\mathrm{sel}}(\omega )P\omega +{\displaystyle \frac{|P\omega (X)|}{\omega }}\omega `$
$``$ $`2P\omega \mathrm{\hspace{0.33em}2}P\omega `$
where $`P:=sup_{\omega 1}P\omega <\mathrm{}`$ because $`P`$ is a bounded operator by assumption, and clearly $`|P\omega (X)||P\omega |(X)=P\omega `$.

Let now $`\omega ,\omega ^{}`$ both be non-zero. Then

$$\mathrm{\Phi }_{\mathrm{sel}}(\omega )\mathrm{\Phi }_{\mathrm{sel}}(\omega ^{})P\omega \omega ^{}+\frac{P\omega ^{}(X)}{\omega ^{}}\omega ^{}\frac{P\omega (X)}{\omega }\omega .$$

Observe that $`\frac{P\omega (X)}{\omega }=P\left(\frac{\omega }{\omega }\right)(X)`$. The second term on the right hand side of the above equation is then clearly majorized by $`P\omega ^{}\omega +c\omega `$ where

$`c`$ $`=`$ $`\left|P\left({\displaystyle \frac{\omega ^{}}{\omega ^{}}}\right)(X)P\left({\displaystyle \frac{\omega }{\omega }}\right)(X)\right|=\left|P\left({\displaystyle \frac{\omega ^{}}{\omega ^{}}}{\displaystyle \frac{\omega }{\omega }}\right)(X)\right|`$
$``$ $`P\left({\displaystyle \frac{\omega ^{}}{\omega ^{}}}{\displaystyle \frac{\omega }{\omega }}\right){\displaystyle \frac{P}{\omega \omega ^{}}}\omega \omega ^{}\omega ^{}\omega .`$
Next, observe that

$$\omega \omega ^{}\omega ^{}\omega \omega ^{}\left|\omega \omega ^{}\right|+\omega ^{}\omega ^{}\omega \mathrm{\hspace{0.33em}2}\omega ^{}\omega \omega ^{}$$

so that we finally get

$$\mathrm{\Phi }_{\mathrm{sel}}(\omega )\mathrm{\Phi }_{\mathrm{sel}}(\omega ^{})\mathrm{\hspace{0.33em}4}P\omega \omega ^{}.$$

Together with the previous calculation, we see that $`\mathrm{\Phi }_{\mathrm{sel}}`$ is globally Lipschitz, with Lipschitz constant $`4`$. $`\mathrm{}`$

So, we know that the IPL equation (30) defines a unique flow. As before, we have to check what happens with $`_+(X)`$ under the semiflow in forward time. Since $`\omega _0=0`$ trivially implies $`\omega _t0`$ for all $`t0`$, we exclude this case from now on. Note that $`\omega _00`$ results in $`\omega _t>0`$ for all $`t0`$, due to uniqueness. Let $`\omega _+(X)`$ be a positive measure and $`E`$ a Borel set such that $`\omega (E)=0`$. This implies $`\mathrm{\Phi }_{\mathrm{sel}}(\omega )(E)=(P\omega )(E)0`$ because $`P`$ itself satisfies the positive minimum principle (A2) by assumption. Also, for any $`\omega _+(X)`$, we have

$$\mathrm{\Phi }_{\mathrm{sel}}(\omega )(X)=P\omega (X)\frac{P\omega (X)}{\omega }\omega (X)=\mathrm{\hspace{0.33em}0}$$

because $`\omega (X)=\omega `$ for positive measures. Together with Fact 6, we see that assumptions (A1) – (A3) are satisfied, and we can invoke Theorem 1.

###### Proposition 9

Assume that the linear operator $`P`$ is bounded and satisfies (A2). Then the abstract Cauchy problem of the IPL equation $`(\text{30})`$ with initial condition $`\omega _0`$ has a unique solution. The cone of positive measures is positive invariant under the flow, and the norm of positive measures is preserved. In particular, $`𝒫(X)`$ is positive invariant. $`\mathrm{}`$

Remark: We would like to mention that the assumption of bounded $`P`$ is somewhat restricted. For non-compact $`X`$, many interesting selection models lead to unbounded $`P`$. For mutation and selection alone, the more general situation has been investigated in  and, more recently, in , in the framework of analytic semigroups, compare \[19, Ch. II.4.a\]. Our emphasis here is on the basic structure that emerges from the interaction with recombination; this will also carry over to more general cases.

Before we proceed, let us make the following observation.

###### Fact 7

If the linear operator $`P`$ is bounded and satisfies the positive minimum principle, the same is true of $`P^{}=P+c\mathbf{\hspace{0.17em}1}`$ for arbitrary $`c`$. Furthermore, the flow of the IPL equation $`(\text{30})`$ on $`_+(X)`$ remains unchanged if $`P`$ is replaced by $`P^{}`$.

Proof: If $`\nu `$ is a positive measure and $`E`$ a Borel set with $`\nu (E)=0`$, then $`P^{}\nu (E)=P\nu (E)+c\nu (E)=P\nu (E)0`$ because $`P`$ satisfies (A2) by assumption. Since $`P^{}`$ is still bounded, the IPL equation (30) with $`P^{}`$ in place of $`P`$ conforms to Proposition 9. If $`\omega _+(X)`$, we obtain

$$P^{}\omega \frac{P^{}\omega (X)}{\omega }\omega =P\omega +c\omega \frac{P\omega (X)}{\omega }\omega \frac{c\omega (X)}{\omega }\omega =P\omega \frac{P\omega (X)}{\omega }\omega $$

from which the claim follows. $`\mathrm{}`$

Once again, although the ODE (30) is nonlinear, it can be solved in closed terms. This time, we employ Thompson’s trick  through the substitution

$$\eta _t=\vartheta (t)\omega _t$$
(32)
where $`\omega _t`$ is a solution of (30). One then obtains

$$\dot{\eta }_t=\left(\dot{\vartheta }(t)\vartheta (t)\frac{P\omega _t(X)}{\omega _t}\right)\omega _t+P\eta _t$$

and a significant simplification is reached if the term in brackets vanishes because the remaining ODE is then linear. This is achieved by the choice

$$\vartheta (t)=\mathrm{exp}\left(_0^t\frac{P\omega _\tau (X)}{\omega _\tau }d\tau \right)=\mathrm{exp}\left(\frac{1}{\omega _0}_0^tP\omega _\tau (X)d\tau \right)$$
(33)
where the second step follows from Proposition 9. Clearly, $`\vartheta (t)`$ is well defined (whenever $`\omega _00`$, which is all we need), and we have reduced the Cauchy problem of (30) to that of the simple linear evolution equation

$$\dot{\eta }=P\eta .$$
(34)
This ODE defines a uniformly continuous positive semigroup (since $`P`$ was assumed to be bounded and to satisfy (A2), the positive minimum principle). The solution of (34) will no longer have fixed norm, but one can always get back to $`\omega _t`$ via

$$\omega _t=\frac{\eta _0}{\eta _t}\eta _t.$$

Note that $`\eta _0=\omega _0`$ and $`\omega _t\eta _0`$.

Let us next consider the function

$$L(t)=\frac{P\omega _t(X)}{\omega _t}$$
(35)
which is defined on any orbit of the flow of (30). $`L(t)`$ is of particular interest on orbits of positive measures, where it admits the interpretation as mean (or averaged) fitness. Here, we know $`\omega _t\omega _0`$ by Proposition 9, so that we obtain

$`\dot{L}(t)`$ $`=`$ $`{\displaystyle \frac{\frac{\mathrm{d}}{\mathrm{d}t}\left(P\omega _t(X)\right)}{\omega _t}}={\displaystyle \frac{P\dot{\omega }_t(X)}{\omega _t}}={\displaystyle \frac{1}{\omega _t}}P\left(P\omega _t{\displaystyle \frac{P\omega _t(X)}{\omega _t}}\omega _t\right)(X)`$
$`=`$ $`P^2\left({\displaystyle \frac{\omega _t}{\omega _t}}\right)(X)\left(P\left({\displaystyle \frac{\omega _t}{\omega _t}}\right)(X)\right)^2`$
which has the form of a variance. So we can state

###### Proposition 10

If, under the assumptions of Proposition 9, $`P`$ satisfies the condition $`(P\omega (X))^2P^2\omega (X)`$ on $`𝒫(X)`$, the function $`L(t)`$ of $`(\text{35})`$ is a Lyapunov function for the flow of $`(\text{30})`$ on the positive cone $`_+(X)`$.

Proof: From the above calculation, it is clear that $`\dot{L}(t)0`$ on all orbits in $`_+(X)`$ if $`P`$ satisfies the inequality $`(P\omega (X))^2P^2\omega (X)`$ on $`𝒫(X)`$, so $`L`$ cannot decrease along such an orbit. $`\mathrm{}`$

Remark: Our definition of a Lyapunov function on $`_+(X)`$ is global and (up to a sign) that of \[4, Ch. 18\]. Note that the stricter version of , where $`\dot{L}(t)=0`$ would correspond to a unique equilibrium on $`_+(X)`$, is not so useful here because the asymptotic state (as $`t\mathrm{}`$) of the selection equation depends on the initial condition, i.e. there is no unique equilibrium in general. However, one might profit from the use of local Lyapunov functions, compare \[26, Thm. 1.0.2 (iii)\], but we do not expand on this here.

The condition on $`P`$ can be reformulated by noting that

$$P^2\omega (X)(P\omega (X))^2=(Pc\mathbf{\hspace{0.17em}1})^2\omega (X)$$

with $`c=P\omega (X)`$. A sufficient condition for Proposition 10 to hold is then that $`(Pc\mathbf{\hspace{0.17em}1})^2`$ is a positive operator for all (or sufficiently many) $`c`$. A particularly well studied case of this is when $`X`$ is finite and $`P`$ is a diagonal matrix in the canonical basis consisting of the extremal measures of $`𝒫(X)`$. In this case, Proposition 10 is known as Fisher’s fundamental theorem, see  for details. In the more general case, Lyapunov functions may be considered even more important since they determine the ‘direction’ of the evolution process in a situation where little information is available otherwise, since the solution given by Eq. (36) is not very explicit then.

The results of this section can also be formulated for the (sub-)space $`^{}`$ of $`(X)`$, if it is invariant under the action of $`P`$. In view of the product structure of $`X`$, let us now assume that we have $`P=_{i=0}^nP_i`$, with bounded $`P_i`$ that are locally represented by $`p_i`$ (as with $`q_i`$ versus $`Q_i`$ before). Clearly, $`P`$ maps $`^{}`$ into itself. We call this situation additivity across sites, in complete analogy to our previous discussion of mutation. We can then rewrite our solution as

$$\eta _t=\mathrm{exp}(tP)\eta _0=\left(\underset{i=0}{\overset{n}{}}\mathrm{exp}(tp_i)\right)\eta _0.$$
(36)
With some further restrictions on the linear operator $`P`$, an analogue of Proposition 10 remains true even in the presence of recombination. This rests on the applicability of Lemma 5. We thus consider the IPL equation

$$\dot{\omega }=P\omega \frac{P\omega (X)}{\omega }\omega +\underset{\alpha L}{}\varrho _\alpha (R_\alpha \mathrm{𝟏})(\omega )=\mathrm{\Phi }_{\mathrm{sel}}(\omega )+\mathrm{\Phi }_{\mathrm{rec}}(\omega )$$
(37)
whose Cauchy problem has all the nice properties we need, see Proposition 11 below in the special case $`Q=0`$. We now assume:

1. $`P`$ has complete product structure as a generator, i.e. $`P=_{i=0}^nP_i`$ with $`P_i=\mathrm{𝟏}\mathrm{}\mathrm{𝟏}p_i\mathrm{𝟏}\mathrm{}\mathrm{𝟏}`$.
2. Each $`P_i`$ is itself a bounded, strictly positive operator.

If $`\omega _t`$ is a solution of (37), we again define $`L(t)`$ as in (35) and obtain, by Lemma 5,

$`\dot{L}(t)`$ $`=`$ $`{\displaystyle \frac{1}{\omega _t}}P\left(P\omega _t{\displaystyle \frac{P\omega _t(X)}{\omega _t}}\omega _t+{\displaystyle \underset{\alpha L}{}}\varrho _\alpha (R_\alpha \mathrm{𝟏})(\omega _t)\right)(X)`$
$`=`$ $`\left[P^2\left({\displaystyle \frac{\omega _t}{\omega _t}}\right)(X)\left(P\left({\displaystyle \frac{\omega _t}{\omega _t}}\right)(X)\right)^2\right]+{\displaystyle \frac{_{i=0}^n\mathrm{\Phi }_{\mathrm{rec}}(P_i\omega _t)(X)}{\omega _t}}.`$
The last term vanishes due to our general assumptions because $`P_i\omega _t>0`$ and then $`\mathrm{\Phi }_{\mathrm{rec}}(P_i\omega _t)(X)=0`$ due to (A3). So, we are back to the condition already encountered above. To summarize:

###### Theorem 5

Let $`P=_{i=0}^nP_i`$ satisfy the assumptions of Proposition 9, and let each $`P_i`$ be a bounded, strictly positive operator with complete product structure. If $`P`$ also satisfies the condition $`(P\omega (X))^2P^2\omega (X)`$ on $`𝒫^{}`$, the function $`L(t)`$ of $`(\text{35})`$ is a Lyapunov function for the flow of $`(\text{37})`$ on the positive cone $`_+^{}`$. $`\mathrm{}`$

In the absence of recombination, there are other Lyapunov functions known for certain combinations of selection with mutation. They rely on the spectral theorem applied to $`P+Q`$, see . In selection-recombination equations where $`P`$ violates the product structure, the mean fitness $`L(t)`$ need no longer be a Lyapunov function. Moreover, the possibility of periodic solutions  demonstrates that, in more general (diploid) models (e.g. with dominance), no meaningful Lyapunov function is to be expected.

## 6 All three

In this last step, we combine all three processes, with the general assumptions as before. In view of the inherent product structure, we only consider the dynamics on the Banach space $`^{}`$. The IPL equation now reads

$`\dot{\omega }`$ $`=`$ $`\mathrm{\Phi }_{\mathrm{mut}}(\omega )+\mathrm{\Phi }_{\mathrm{rec}}(\omega )+\mathrm{\Phi }_{\mathrm{sel}}(\omega )`$
$`=`$ $`(Q+P)\omega {\displaystyle \frac{P\omega (X)}{\omega }}\omega +{\displaystyle \underset{\alpha L}{}}\varrho _\alpha \left(R_\alpha \mathrm{𝟏}\right)(\omega )`$
and we immediately get the following result, again from Theorem 1, and Lemma 4.

###### Proposition 11

Let $`Q`$ be a bounded Markov generator and $`P`$ a bounded generator of a positive semigroup, both of product form. Let $`R_\alpha `$ be the recombinators of Eq. $`(\text{9})`$. Then, the abstract Cauchy problem of the IPL equation $`(\text{6})`$ has a unique solution. The cone $`_+^{}`$ is positive invariant and the norm of positive measures is preserved under the forward flow. $`\mathrm{}`$

Remark: Since $`Q`$ is a Markov generator, we know from Section 2 that $`Q\omega (X)=0`$ for all $`\omega _+^{}`$, and we could also start from an IPL equation where $`Q`$ is absorbed into $`P`$ — it would give the same flow on $`_+^{}`$. We retain the separation into mutation and selection because, in more general situations, it is often adequate from both the biological and the mathematical point of view (for example, the mutation operator is usually bounded, but the selection operator may be unbounded); for review, see \[11, Ch. IV\]. We will also combine $`Q`$ and $`P`$, but only after Thompson’s linearization transformation.

Let $`\omega _t`$, $`t0`$, be the solution for initial condition $`\omega _0`$. Define $`\eta _t`$ as above in (32), with $`\vartheta (t)`$ of (33). Then, $`\omega _t`$ is a solution of (6) if and only if $`\eta _t`$ solves the reduced IPL equation

$$\dot{\eta }=S\eta +\mathrm{\Phi }_{\mathrm{rec}}(\eta )$$
(39)
where $`S=Q+P`$ is the bounded generator of a uniformly continuous semigroup of positive operators. Note that the right hand side of (39) still satisfies assumptions (A1) and (A2), but no longer (A3). So, the corresponding Cauchy problem still has a unique solution, with $`_+^{}`$ being positive invariant, but the norm of positive measures need no longer be preserved under the flow in forward time — and this is precisely the point of this exercise!

From now on, we generally assume that both mutation and selection are adapted to the special product form of our state space $`X`$, so $`S=_{i=0}^nS_i`$ (with corresponding local operator $`s_i`$). Hence, $`\mathrm{exp}(tS)`$ is again a tensor product of local operators.

###### Lemma 6

If $`S=_{i=0}^nS_i`$ is the bounded generator of a uniformly continuous semigroup of positive operators, then we have $`\mathrm{exp}(tS)R_\alpha =R_\alpha \mathrm{exp}(tS)`$ on $`_+^{}`$, for all $`t0`$ and $`\alpha L`$.

Proof: Fix $`t0`$ and set $`W=\mathrm{exp}(tS)`$. This is a positive operator by assumption. Also, since $`S`$ is bounded, $`\nu >0`$ implies $`\mathrm{exp}(tS)\nu >0`$ and $`W`$ is strictly positive. The result then follows from Lemma 5. $`\mathrm{}`$

This result means that we can use all our above methods again and construct immediately the solution of (39). At this point, we particularly profit from our approach in that we can still solve the case with (additive) selection. In the context of Haldane linearization, any form of selection has, so far, appeared as a major obstacle, due to the fact that the flow induced by $`P`$ fails to preserve the norm of positive measures .

###### Theorem 6

If $`S=_{i=0}^nS_i`$ satisfies the assumptions of Lemma 6, the solution of the reduced IPL equation $`(\text{39})`$, with initial condition $`\eta _0_+^{}`$, is given by

$$\eta _t=\mathrm{exp}(tS)\underset{GL}{}a_G(t)R_G(\eta _0)$$

with the coefficients $`a_G(t)`$ of $`(\text{13})`$. The solution of the abstract Cauchy problem for the original IPL equation $`(\text{6})`$ emerges from here via

$$\omega _t=\frac{\omega _0}{\eta _t}\eta _t$$

where $`\omega _0=\eta _0`$. If $`\omega _0𝒫(X)`$, then $`\{\omega _tt0\}`$ is a one-parameter family of probability measures. $`\mathrm{}`$

In line with our previous reasoning, we can determine the asymptotic behaviour,

$$\eta _t\underset{i=0}{\overset{n}{}}\left(\mathrm{exp}(ts_i)(\pi _i.\eta _0)\right),$$

where we have used the product structure of $`\mathrm{exp}(tS)`$ and the fact that all neglected terms, as $`t\mathrm{}`$, are exponentially small in comparison. The meaning for $`\omega _t`$ is, once again, that stationary measures are complete product measures, and the properties of the linear operators $`s_i`$ determine whether there is a unique global equilibrium measure. This is connected to the general Perron-Frobenius theory of positive operators which is rather involved in general, see \[42, Ch. V.5\] and . If, however, $`X`$ is finite (so that $`^{}`$ is finite-dimensional, and $`^{}=(X)`$) and all $`s_i`$ are irreducible, there are unique $`\nu _i𝒫(X_i)`$ so that $`\mathrm{exp}(ts_i)\nu _i=\mathrm{exp}(t\lambda _i)\nu _i`$ with $`\lambda _i`$ being the largest eigenvalue of $`s_i`$. In this case, as a simple calculation shows, we obtain

$$\omega _t\nu _0\mathrm{}\nu _n$$

in the $`.`$-topology, as $`t\mathrm{}`$, for any initial condition $`\omega _0𝒫(X)`$.

Also, the following observation results immediately from Theorem 6.

###### Corollary 3

If an initial condition $`\omega _0𝒫^{}`$ is a product measure at link $`\alpha L`$, this is also true of the corresponding solution $`\omega _t`$ of $`(\text{6})`$, for all $`t0`$. In particular, if $`\omega _0`$ is a complete product measure, this remains the case under the forward flow, i.e. for all $`\omega _t`$ with $`t0`$. $`\mathrm{}`$

Let us return to the general discussion. The remainder is then a copy of what we did in Section 4, with $`Q`$ replaced by $`S`$. In particular, we get

$$\eta _t=\underset{KL}{}b_K(t)T_K\left(\mathrm{exp}(tS)\eta _0\right)$$

from which one can, once again, determine the linkage disequilibria. Note, however, that the meaning has changed now, because the norm of $`\eta _t`$ varies with time. In particular, one has to consider the quotient $`\eta _t/\eta _t`$, rather than $`\eta _t`$ alone, to extract the correct behaviour for the linkage disequilibria $`F_G^t(j_1,\mathrm{},j_k)=T_G(\omega _t)\left(j_1,\mathrm{},j_k\right)`$. To be concrete, observe first that

$$\frac{d}{dt}T_G(\eta _t)=\left(S\underset{\alpha \overline{G}}{}\varrho _\alpha \right)T_G(\eta _t)$$

in perfect analogy with (29). Since $`T_G`$ is positive homogeneous of degree one (Eq. (22)), and $`\eta _t=\eta _t(X)`$ for positive measures, one obtains

$$\frac{d}{dt}T_G(\omega _t)=\left(S\frac{S\omega _t(X)}{\omega _t}\mathbf{\hspace{0.17em}1}\underset{\alpha \overline{G}}{}\varrho _\alpha \right)T_G(\omega _t).$$

Clearly, knowledge of the mean fitness, $`S\omega _t(X)/\omega _t`$, is now required to determine the dynamics of the linkage disequilibria.

## 7 Afterthoughts

In this article, we have constructed an explicit solution of the single-crossover recombination model in continuous time, with mutation and additive selection. It is quite astonishing that such a solution should be possible at all – after all, it is an explicit representation of a nonlinear semigroup. However, it is no coincidence that this works in continuous time, rather than in discrete time. Let us discuss this for recombination alone. The discrete-time analogue of our single-crossover model is the so-called model with complete interference :

$$\omega _{n+1}=\underset{\alpha L}{}\varrho _\alpha R_\alpha (\omega _n)+\left(1\underset{\alpha L}{}\varrho _\alpha \right)\omega _n.$$
(40)
Similar as this may look to its continuous-time relative, the probabilistic structure is quite different. Single crossovers in continuous time imply independence of links, as expressed in the coefficient functions (13) and the resulting factorization property (Lemma 2). In contrast, a second crossover is inhibited for the duration of an entire generation in discrete time, due to interference of crossovers with each other (hence the name); see also . As a result, independence is lost, which makes the discrete model inherently more difficult.

Of course, this also applies to the situation with selection. Models of recombination and selection based on independent sites and finite site spaces have been thoroughly investigated in the population genetics literature, see  for some key references and  for recent comprehensive reviews. Independence of sites with respect to selection is reflected by a tensor product structure of $`P`$, may be interpreted as lack of interaction between genes, and is known as absence of epistasis in genetics. More precisely, since the dynamical systems mostly considered so far were in discrete time, a comparison with our setting is more adequate at the level of the semigroup, rather than that of the generator.

Two notions of independence have been used, compare , which would translate into our setting as either $`\mathrm{exp}(P)=_i\mathrm{exp}(P_i)=_i\mathrm{exp}(p_i)`$ (‘multiplicative fitness’) or as $`\mathrm{exp}(P)`$ replaced by $`_i\mathrm{exp}(P_i)`$ (‘additive fitness’). Previously, much emphasis has been on the effects of dominance (i.e. the interaction between the two alleles joined in a diploid genotype). This may lead to multiple equilibria, which need not all be of product type, and astonishing differences in the qualitative behaviour of the multiplicative and additive scenario are observed, as reviewed in . However, these effects are absent if there is no dominance (as in our model); in particular, all equilibria are then of product type. Thus, our simple continuous-time model might well serve as an exactly solved reference case which also captures the qualitative features of the corresponding models in discrete time, although no explicit solution is available there.

Now, the logical next step would be to extend the analysis to the inclusion of interactions between sites, which occur as soon as selection is no longer additive across sites. Alas, this is much more involved, and even the simplest cases go far beyond what we have outlined above. The reason is that selection now forces the introduction of further terms in the right hand side of the IPL equation so that the corresponding semigroups no longer commute with recombination. Nevertheless, several situations can be envisioned that admit at least a perturbative approach. In line with the single-crossover assumption, an expansion for small recombination rates would be appropriate, in contrast to the well-known quasi-linkage-equilibrium approach for large recombination rates (for review, see ). We hope to report on some progress in this direction soon.

### Appendix: Moments versus correlations

As mentioned above, it is often desirable to separate effects that stem from mutual interactions of differently many “particles” or, as in the above discussion, from specification at a different number of sites. For two sites, correlation $`C`$ and moments $`F`$ are related by $`C(\{i,j\})=F(\{i,j\})F(\{i\})F(\{j\})`$, where the arguments are meant as symbolic labels. Since this is a rather general structure, we briefly describe its systematic treatment by means of Möbius inversion, also known as inclusion-exclusion principle.

Let $`S=\{1,2,\mathrm{},k\}`$ be a finite set which will serve as the index set of the particles or the specified sites, the latter through $`j_1,\mathrm{},j_{\mathrm{}}`$. Let $`𝒜=\{A_1,\mathrm{},A_p\}`$ be a partition of $`S`$, i.e. $`S`$ is the disjoint union of the members of $`𝒜`$. Unlike before, the partition need not be ordered. Let the partition $`=\{B_1,\mathrm{},B_q\}`$ be a refinement of $`𝒜`$, so that

$$A_1=B_{j_{1,1}}\mathrm{}B_{j_{1,n_1}},\mathrm{},A_p=B_{j_{p,1}}\mathrm{}B_{j_{p,n_p}}$$

where $`\{\{j_{1,1},\mathrm{},j_{1,n_1}\},\mathrm{},\{j_{p,1},\mathrm{},j_{p,n_p}\}\}`$ is a partition of $`\{1,\mathrm{},q\}`$, hence $`n_1+\mathrm{}+n_p=q`$. We write $`𝒜`$ in this case, where $``$ defines a partial order which makes $`S`$ into a poset. The corresponding Möbius function, compare \[10, p. 86\], is given by

$`\mu (,𝒜)`$ $`=`$ $`{\displaystyle \underset{i=1}{\overset{p}{}}}(1)^{n_i1}(n_i1)!`$
$`=`$ $`(1)^{p+n_1+\mathrm{}+n_p}(n_11)!\mathrm{}(n_p1)!`$
If $`𝒞`$ is any refinement of $`𝒜`$, $`\mu `$ satisfies the formula

$$\underset{𝒞\text{\$̣mathcal{B}\$}𝒜}{}\mu (,𝒜)=\{\begin{array}{cc}1\hfill & \text{if }𝒜=𝒞\hfill \\ 0\hfill & \text{otherwise.}\hfill \end{array}$$

Let us now, for a partition $`𝒜=\{A_1,\mathrm{},A_p\}`$, introduce the function $`F(𝒜)=F(A_1)\mathrm{}F(A_p)`$, and similarly for the correlations, $`C`$. These quantities are related by

$$F(𝒜):=\underset{𝒜}{}C()=\underset{𝒜}{}\underset{B}{}C(B)$$
(42)
because this precisely reflects the idea to separate off contributions from subsets of different cardinality. The Möbius inversion formula then gives the following formula for the special case that $`𝒜=\{A\}`$:

$$C(A)=\underset{𝒜}{}F()\mu (,𝒜)=\underset{𝒜}{}(1)^{||1}\left(||1\right)!\underset{i=1}{\overset{||}{}}F(B_i)$$
(43)
where $`||`$ denotes the number of sets in the partition $`=\{B_1,\mathrm{},B_{||}\}`$. The following example might illustrate this formula:

$`C(\{1,2,3\})`$ $`=`$ $`F(\{1,2,3\})+\mathrm{\hspace{0.33em}2}F(\left\{1\right\})F(\left\{2\right\})F(\left\{3\right\})`$
$`F(\left\{1\right\})F(\{2,3\})F(\left\{2\right\})F(\{1,3\})F(\left\{3\right\})F(\{1,2\})`$
which is to be compared with

$`F(\{1,2,3\})`$ $`=`$ $`C(\{1,2,3\})+C(\left\{1\right\})C(\left\{2\right\})C(\left\{3\right\})`$
$`+C(\left\{1\right\})C(\{2,3\})+C(\left\{2\right\})C(\{1,3\})+C(\left\{3\right\})C(\{1,2\})`$
according to (42). Let us finally remark that formula (43) can be applied factorwise if $`𝒜=\{A_1,\mathrm{},A_p\}`$ because then $`C(𝒜)=C(A_1)\mathrm{}C(A_p)`$ by definition.

### Acknowledgements

It is our pleasure to thank Reinhard Bürger, Hans-Otto Georgii, Joachim Hermisson, Achim Klenke and Manfred Wolff for a number of clarifying discussions, and Ulrich Hermisson and Oliver Redner for carefully reading the manuscript. M.B. would like to thank Robert V. Moody and the Department of Mathematical Sciences of the University of Alberta (Edmonton, Canada) for hospitality, where part of this work was done.