The model-specific Markov embedding problem for symmetric group-based models.

We study model embeddability, which is a variation of the famous embedding problem in probability theory, when apart from the requirement that the Markov matrix is the matrix exponential of a rate matrix, we additionally ask that the rate matrix follows the model structure. We provide a characterisation of model embeddable Markov matrices corresponding to symmetric group-based phylogenetic models. In particular, we provide necessary and sufficient conditions in terms of the eigenvalues of symmetric group-based matrices. To showcase our main result on model embeddability, we provide an application to hachimoji models, which are eight-state models for synthetic DNA. Moreover, our main result on model embeddability enables us to compute the volume of the set of model embeddable Markov matrices relative to the volume of other relevant sets of Markov matrices within the model.

Introduction

The embedding problem for stochastic matrices, also known as Markov matrices, deals with the question of deciding whether a stochastic matrix M is the matrix exponential of a rate matrix Q. A rate matrix, also known as a Markov generator, has non-negative non-diagonal entries and row sums equal to zero. If a stochastic matrix satisfies such a property and can be expressed as a matrix exponential of a rate matrix, namely , then M is called embeddable. Applications of the embeddability property vary from biology and nucleotide substitution models (Verbyla et al. 2013) to mathematical finance (Israel et al. 2001). For a first formulation of the embedding problem see Elfving (1937). An account of embeddable Markov matrices is provided in Davies (2010). The embedding problem for matrices is due to Kendall and first published by Kingman (1962), for matrices is fully settled in a series of papers Carette (1995), Chen and Chen (2011), Cuthbert (1973), Israel et al. (2001) and Johansen (1974), while for stochastic matrices has been recently solved in Casanellas et al. (2020b). In general, when the size n of the stochastic matrix is greater than 4, the work Casanellas et al. (2020b) establishes a criterion for deciding whether a generic Markov matrix with distinct eigenvalues is embeddable and proposes an algorithm that lists all its Markov generators. In the present paper we study a refinement of the classical embedding problem, called the model embedding problem for a class of stochastic matrices coming from phylogenetic models.

Phylogenetics is the field that aims at reconstructing the history of evolution of species. A phylogenetic model is a mathematical model used to understand the evolutionary process given genetic data sets. The most popular phylogenetic models are nucleotide substitution models which use aligned DNA sequence data to study the molecular evolution of DNA. A comprehensive treatment of phylogenetic methods is given by Felsenstein, who is considered the initiator of statistical phylogenetics, in his seminal book Felsenstein (2003). Algebraic and geometric methods have been employed with great success in the study of phylogenetic models leading to an explosion of related research work and the establishment of the field of phylogenetic algebraic geometry, also known as algebraic phylogenetics; see Allman and Rhodes (2003), Baños et al. (2016), Casanellas and Fernández-Sánchez (2007), Cavender and Felsenstein (1987), Gross and Long (2018), Evans and Speed (1993), Lake (1987), Pachter and Sturmfels (2004) and Sturmfels and Sullivant (2005) for a non-exhaustive list of publications.

To build such a phylogenetic model, we first require a phylogenetic tree T, which is a directed acyclic graph comprising of vertices and edges representing the evolutionary relationships of a group of species. The vertices with valency 1 are called the leaves of the tree. The tree is considered rooted and the direction of evolution is from the root towards the leaves. On each vertex of the tree T, we associate a random variable with k possible states, where in phylogenetics k is often taken to be 2, for the binary states , or 4, to represent the four types of DNA nucleotides . We also require a transition matrix (also known as a mutation matrix) corresponding to each edge e of the tree, where the entries of this matrix represent the probabilities of transition between states. In a phylogenetic tree, the leaves correspond to extant species and so the random variables at the leaves are observed, while the interior vertices correspond to possibly extinct species and so the random variables at the interior vertices are hidden.

In this paper we are focusing on symmetric group-based substitution models. Substitution models are a class of phylogenetic models which use a Markov process to describe the substitution of nucleotides over time in a given DNA sequence and for which the transition matrices along an edge e are stochastic matrices of the form . Group-based models are a special class of substitution models, in which the matrices can be pairwise distinct, but they can all be simultaneously diagonalizable by a linear change of coordinates given by the discrete Fourier transform of an abelian group, also called a commutative group. For example, the Cavender–Farris–Neyman (CFN) model (Cavender 1978; Farris 1973; Neyman 1971), as well as the Jukes–Cantor (JC) (Jukes and Cantor 1969), the Kimura-2 parameter (K2P) (Kimura 1980) and the Kimura-3 parameter (K3P) (Kimura 1981) models for DNA are all group-based phylogenetic models. In Sturmfels and Sullivant (2005), it is established that through the discrete Fourier transform group-based models correspond to toric varieties, which are geometric objects with nice combinatorial properties. We are interested in symmetric group-based substitution models. Namely, apart from distinct and simultaneously diagonalizable, the transition matrices are also symmetric square matrices. The symmetricity assumption guarantees that the eigenvalues of rate and transition matrices of a group-based model are real, a property that we use in the proof of our main theorem. Symmetric models are a subset of a special class of models called time-reversible models, where the Markov process appears identical when moving forward or backward in time. Our results apply to group-based models following an ergodic time-reversible Markov process, as in this case the rate matrices Q are symmetric according to Pachter and Sturmfels (2004, Lemma 17.2).

The classical embedding problem is concerned with deciding which square Markov matrices are embeddable, namely given a Markov matrix M whether there exists a rate matrix Q such that . A variant of the embedding problem that asks for a reversible Markov generator for a stochastic matrix is studied in Jia (2016). When we impose the assumption that the rate matrix Q follows the corresponding model conditions, we arrive at a different refined notion of embeddability called model embeddability. The embeddability of circulant and equal-input stochastic matrices is studied in Baake and Sumner (2020). In the current paper, we focus on the -embeddability for matrices corresponding to symmetric group-based substitution models. The -embeddability means that we require that the rate matrices Q preserve the symmetric group-based structure imposed by the abelian group and the symmetric labelling L, which we define at the beginning of Sect. 2. Model embeddability for symmetric group-based models is relevant both for homogeneous and inhomogeneous time-continuous processes, as group-based models are Lie Markov models, and hence multiplicatively closed (Sumner et al. 2012; Verbyla et al. 2013). A study of the set of embeddable and model-embeddable matrices corresponding to the Jukes–Cantor, Kimura-2 and Kimura-3 DNA substitution models, which are all symmetric group-based models, is undertaken in Casanellas et al. (2020a) and Roca-Lacostena and Fernández-Sánchez (2018). In particular, a full characterisation of the set of embeddable Kimura 2-parameter matrices is provided in Casanellas et al. (2020a), which together with the results of Roca-Lacostena and Fernández-Sánchez (2018) fully solve the embedding problem for the Kimura 3-parameter model as well. Although model embeddability, which is a refined notion of embeddability imposed by the model structure, implies classical embeddability, the converse is generally not true (see also Roca-Lacostena and Fernández-Sánchez 2018, Example 3.1).

The main result of this paper is a characterization of -embeddability for any abelian group equipped with a symmetric -labeling function L in Theorem 1. We provide necessary and sufficient conditions which the eigenvalues of the stochastic matrix of the model need to satisfy for the matrix to be -embeddable. To showcase our result, we first introduce three group-based models with the underlying group , based on the hachimoji DNA system introduced in Hoshika et al. (2019). Hachimoji, a Japanese word meaning “eight letters”, is used to describe a synthetic analog of the nucleic acid DNA, where we have the four natural nucleobases and furthermore an additional four synthetic nucleotides . We then apply Theorem 1 to characterise the model embeddability for the three hachimoji DNA models. The three models are called hachimoji 7-parameter, hachimoji 3-parameter and hachimoji 1-parameter models, which can be thought of as generalisations of the Kimura 3, Kimura 2 and Jukes–Cantor models respectively. Finally, the characterisation of model embeddability in terms of eigenvalues enables us to compute the volume of the -embeddable Markov matrices and compare this volume with volumes of other relevant sets of Markov matrices. For the general Jukes–Cantor model, which includes the hachimoji 1-parameter model, the volumes can be derived exactly; for the hachimoji 3-parameter model and for the hachimoji 7-parameter model symbolically and numerically.

The outline of the paper is the following. Section 2 gives a mathematical background covering notions such as the labeling functions, group-based models and the discrete Fourier transform. Section 3 introduces symmetric -compatible labelings which is a class of labeling functions with particularly nice properties. Section 4 presents the main result of this paper about the model embedding problem for symmetric group-based models equipped with a certain labeling function. Then in Sect. 5 we focus on the hachimoji DNA and provide exact characterisation of the model embeddability in terms of eigenvalues of the Markov matrix for the hachimoji 7-parameter, the hachimoji 3-parameter and the hachimoji 1-parameter models. Finally, Sect. 6 presents results on the volume of stochastic matrices that are -embeddable for the three hachimoji group-based models. The code for the computations in this paper is available at https://github.com/ardiyam1/Model-Embeddability-for-Symmetric-Group-Based-Models.

Preliminaries

In this section, we give background on group-based models and the discrete Fourier transform.

Definition 1

Let be a finite additive abelian group and a finite set. A labeling function is any function .

In the group-based model with underlying finite abelian group , states are in bijection with the elements of the group . Fundamental in the definition of a group-based model associated to a finite additive abelian group and a labeling function L is that the rate of mutation from a state g to state h depends only on : That is, the entries of a rate matrix Q arefor a vector satisfying , for all non-zero and , whenever . We say that such Q is a -rate matrix. In this paper, the rate matrices in group-based models are assumed to be symmetric, i.e., for every . Since the matrix exponential of a symmetric matrix is again symmetric, then the entries of a transition matrix arefor a vector satisfying , for all and for all . As we see in Example 2, in general it is not true that whenever . In Sect. 3, we introduce -compatible labeling functions which guarantee this property and then we say that P is a -Markov matrix.

Example 1

Let and . We identify nucleotides with the group elements of as . The Kimura 3-parameter, the Kimura 2-parameter and the Jukes–Cantor models correspond to the following labeling functionsrespectively. The Kimura 3-parameter rate and transition matrices have the formIn the case of the Kimura 2-parameter model additionally , after choosing the ordering of nucleotide bases. This will be further explained in Example 3. The Jukes–Cantor model is the submodel when .

Example 2

Let and L be a labeling function such that and . Consider the -rate matrixIn this rate matrix, and . By direct computation, we getThe matrix P is not a -Markov matrix, since . The entries of P induce a labeling function such that and In this case, the matrix P is a -Markov matrix.

Example 2 shows that the matrix exponential does not necessarily preserve the labeling function associated to a rate matrix. Conversely, Example 3.1 of Roca-Lacostena and Fernández-Sánchez (2018) suggests that a Kimura 3-parameter Markov matrix can be embeddable, despite the fact that it does not have any Markov generator satisfying Kimura 3-parameter model constraints.

Let denote the multiplicative group of complex numbers without zero. A group homomorphism from to is called a character of . The characters of form a group under multiplication, called the character group of and denoted by . Here the product of two characters of the group is defined by for every . The character group is isomorphic to . Given a group isomorphism between and , we will denote by the image of . For a finite group , the values of characters are roots of unity.

Lemma 1

(Pachter and Sturmfels 2005, Lemma 17.1) Let and . Then and , where denotes the complex conjugate of .

Given a function , its discrete Fourier transform is a function defined by

Lemma 2

(Matsen 2008, Section 2) For any real-valued function , the identity holds for all . Moreover, if for all , then for all and is a real-valued function.

In the proof of Theorem 1, we will use that and are real-valued. For this reason, in this paper we consider only group-based models that are equipped with symmetric labeling functions, i.e. for all . In other words, a symmetric group-based model assumes that the transition matrices are real symmetric matrices.

The discrete Fourier transform is a linear endomorphism on . We will denote its matrix by K. In particular, the entries of K are for . The matrix K is symmetric for any finite abelian group (Luong 2009, Section 3.2). The inverse of the discrete Fourier transformation matrix is , where denotes the adjoint of K (Luong 2009, Corollary 3.2.2).

The following lemma describes the relation between functionals f and in the case and for all .

Lemma 3

(Matsen 2008, Lemma 2.2) Let Q be determined by and P be determined by as described earlier in this section such that . Furthermore, assume that and for all . Then, for all .

Lemma 4

Let Q be determined by and P be determined by as described earlier in this section. Furthermore, assume that and for all . Let denote the column of the discrete Fourier transform matrix labeled by g. The eigenpairs of Q (resp. P) are (resp. ) for .

Proof

This result is stated in the proof of  Matsen (2008, Lemma 2.2).

In particular, in the case of a Markov matrix, the column vector of ones is an eigenvector with eigenvalue one. In the case of a rate matrix, the column vector of ones is an eigenvector with eigenvalue zero. A direct consequence of Lemma 4 is that Q and P are diagonalizable by K, i.e. and where and are diagonal matrices with diagonals given by the vectors and of respectively.

-compatible labeling functions

In this section, we introduce a class of labeling functions with the property that the symmetries of the probability vector are preserved under the discrete Fourier transformation. This property is required for any result that is proven using the discrete Fourier transform. Notably, labeling functions for all common group-based models (CFN, K3P, K2P, and JC models) are -compatible.

Definition 2

Let be a finite additive abelian group, a set and a labeling function. Let K be the discrete Fourier transformation matrix for and be the column vector of length whose g-th component is the indeterminate . We say that L is a -compatible labeling function if for every with , we have that and . Here denotes the row of M indexed by group element a.

A labeling function that maps every group element to a different label is trivially -compatible.

Remark 1

In the definition of a -compatible labeling, we require that the matrices K and preserve the symmetries of the vector of labels . For symmetric group-based models, it is enough to require that only K or preserves the symmetries of the vector of labels . Recall thatThe property implies and . If , this means for all . If , then taking into account that gives . Hence .

Remark 2

If a labeling function L is symmetric -compatible and Q is a -rate matrix, then a Markov matrix is a -Markov matrix, i.e. for a vector and whenever .

Example 3

Let . The discrete Fourier transformation matrix for isTo show -compatibility for the three labeling functions from Example 1, it is enough to check that K preserves the symmetries of . The labeling function of the Jukes–Cantor model is -compatible, sinceThe labeling function of the Kimura 2-parameter model is -compatible, sinceIn the literature, usually in the Kimura 2-parameter model. However, here we assume that which is simply due to the fact that we consider the identification . To be more precise, nucleotide bases fall into two categories depending on the molecular mechanisms of the base; purines (A or G) and pyrimidines (C or T). A transition occurs when a purine is substituted by a purine, or a pyrimidine by a pyrimidine. A change from a purine to a pyrimidine, or vice versa, is a transversion. The Kimura 2-parameter model of sequence evolution distinguishes between transitions and transversions to account for the biological fact that transitions occur at higher rate than transversions (Kimura 1980, 1981). The rate and transition matrix of the Kimura 2-parameter model have the formThe reason for choosing this identification and ordering is that we can use the discrete Fourier transform matrix in a format, which better demonstrates that it is the Kronecker product of discrete Fourier transformation matrices for . The labeling function of the Kimura 3-parameter is -compatible, because each group element maps to a different label.

Sturmfels and Sullivant consider a different class of labeling functions, called friendly labelings (Sturmfels and Sullivant 2005). Group-based models with friendly labeling functions are equivalent to -models defined by Michałek (2011, Remark 5.2). -models are constructed using an arbitrary group that has a normal, abelian subgroup which acts transitively and freely on the finite set of states. The importance of -models is that they are toric. We explore connections between friendly labelings and -compatible labelings in “Appendix A”. We conjecture that every symmetric -compatible labeling is a friendly labeling, but not vice versa.

The following lemma provides a necessary condition for -compatible labeling functions.

Lemma 5

Let be a finite additive abelian group, a set and a labeling function. If L is -compatible, then for any .

Let K be the discrete Fourier transformation matrix for . The entries of K are for , which are roots of unity. The row consists of ones. On the other hand, no other row of K consists of ones only, as this would contradict the uniqueness of the identity element in the character group. In particular, every other row of K contains at least one element whose real part is strictly less than one. Thus it is impossible that for .

Table 1 summarizes up to isomorphism all possible symmetric -compatible labeling functions for additive abelian groups of order up to eight. In the table, two group elements receive the same label if they belong to the same subset in a partition of .

Table 1

Symmetric -compatible labelings for abelian groups of order up to isomorphism

n	Group	Symmetric G-compatible labelings
2	Z_2	{{0},{1}}
3	Z_3	{{0},{1,2}}
4	Z_4	{{0},{1,2,3}},{{0},{1,3},{2}}
4	Z_2×Z_2	{{(0,0)},{(0,1),(1,0),(1,1)}},{{(0,0)},{(0,1),(1,0)},{(1,1)}}, {{(0,0)},{(0,1)},{(1,0)},{(1,1)}}
5	Z_5	{{0},{1,2,3,4}}, {{0},{1,4},{2,3}}
6	Z_2×Z_3	{{(0,0)},{(0,1),(0,2),(1,0),(1,1),(1,2)}}, {{(0,0)},{(0,1),(0,2)},{(1,0)},{(1,1),(1,2)}}
7	Z_7	{{0},{1,2,3,4,5,6}},{{0},{1,6},{2,5},{3,4}}
8	Z_8	{{0},{1,2,3,4,5,6,7}},{{0},{1,3,5,7},{2,6},{4}}, {{0},{1,7},{2,6},{3,5},{4}}
8	Z_2×Z_4	{{(0,0)},{(0,1),(0,2),(0,3),(1,0),(1,1),(1,2),(1,3)}}, {{(0,0)},{(0,1),(0,3),(1,1),(1,3)},{(0,2)},{(1,0),(1,2)}}, {{(0,0)},{(0,1),(0,2),(0,3)},{(1,0)},{(1,1),(1,2),(1,3)}}, {{(0,0)},{(0,1),(0,3),(1,0)},{(0,2),(1,1),(1,3)},{(1,2)}}, {{(0,0)},{(0,1),(0,3)},{(0,2)},{(1,0)},{(1,1),(1,3)},{(1,2)}}
8	Z_2×Z_2×Z_2	{{(0,0,0)},{(0,0,1),(0,1,0),(0,1,1),(1,0,0),(1,0,1),(1,1,0),(1,1,1)}}, {{(0,0,0)},{(0,0,1)},{0,1,0),(1,0,0),(1,1,0)},{(0,1,1),(1,0,1),(1,1,1)}}, {{(0,0,0)},{(0,0,1),(1,0,0),(1,0,1)}, {(0,1,0)},{(0,1,1),(1,1,0),(1,1,1)}}, {{(0,0,0)},{(0,0,0),(0,1,0),(0,1,1)},{(1,0,0)},{(1,0,1),(1,1,0),(1,1,1)}}, {{(0,0,0)},{(0,0,1),(0,1,0),(1,0,1),(1,1,0)},{(0,1,1)},{(1,0,0),(1,1,1)}}, {{(0,0,0)},{(0,0,1),(0,1,0),(1,0,0)},{(0,1,1),(1,0,1),(1,1,0)},{(1,1,1)}}, {{(0,0,0)},{(0,0,1),(1,1,1)},{(0,1,0),(0,1,1),(1,0,0),(1,0,1)},{(1,1,0)}}, {{(0,0,0)},{(0,0,1),(0,1,1),(1,0,0),(1,1,0)},{(0,1,0),(1,1,1)},{(1,0,1)}}, {{(0,0,0)},{(0,0,1)},{(0,1,0),(1,0,0)},{(0,1,1),(1,0,1)},{(1,1,0)},{(1,1,1)}}, {{(0,0,0)},{(0,0,1),(0,1,0)},{(0,1,1)},{(1,0,0)},{(1,0,1),(1,1,0)},{(1,1,1}}, {{(0,0,0)},{(0,0,1),(1,0,0)},{(0,1,0)},{(0,1,1),(1,1,0)},{(1,0,1)},{(1,1,1}}, {{(0,0,0)},{(0,0,1)},{(0,1,0)},{(0,1,1)},{(1,0,0)},{(1,0,1)},{(1,1,0)},{(1,1,1)}}

We saw in Example 3 that the labeling function of the Jukes–Cantor model that assigns the same label to each nonzero element of the group is a -compatible labeling. This example can be generalized to other groups.

Lemma 6

Let be a finite abelian group. Let be a labeling function such that and and for . Then the labeling function L is symmetric -compatible.

Clearly the labeling function L is symmetric. Let K be the discrete Fourier transformation matrix for . By Luong (2009, Corollary 3.2.1), we haveHenceHence L is a -compatible labeling function.

We call the model in Lemma 6 the general Jukes–Cantor model. We finish this section with giving another class of labeling functions that are -compatible for every finite abelian group .

Lemma 7

Let be a finite additive abelian group, a finite set and a labeling function such that for any two distinct elements , if and only if . Then L is -compatible.

By Lemma 1, the identity holds for all . ThenThus, the labeling function L is -compatible as is the column vector of indeterminates .

The converse of Lemma 7 is not true in general. Two examples are given by the labeling functions for the Kimura 2-parameter and the Jukes–Cantor model.

Model embeddability

The following theorem is the main result of this paper. It characterizes -embeddable transition matrices in terms of their eigenvalues.

Theorem 1

Fix a finite abelian group , a finite set , and a symmetric -compatible labeling function . Let P be a -Markov matrix. Then P is -embeddable if and only if the vector of eigenvalues of P is in the set

We start by summarizing the idea of the proof. We consider the set that consists of vectors that determine -rate matrices. Our goal is to characterize the set of eigenspectra of Markov matrices that are matrix exponentials of -rate matrices determined by vectors in . The first step is to consider the discrete Fourier transform of the set , which we denote by . By Lemma 4, this set is the set of eigenvalues of the -rate matrices. The second step is to consider the image of the set under coordinatewise exponentiation. This set is precisely , because -rate matrices are diagonalizable by the discrete Fourier transform matrix K by the discussion after Lemma 4 and thus if a -rate matrix Q is determined by thenwhere is the vector of eigenvalues of Q and is the vector of eigenvalues of P.

More specifically, letThe vectors in the set are in one-to-one correspondence with -rate matrices. The image of under the discrete Fourier transform is the setBy Lemma 4, this set is the set of eigenvalues of the -rate matrices.

The image of under the coordinatewise exponentiation is the set of eigenvalues of the -Markov matrices, which we denote by . We claim that is equal to the setIndeed, let . Then because the image of the exponentiation map is positive. The inequality is equivalent to Hence, the equation gives and the inequalities givefor all nonzero . Hence is in the set (4.1). Conversely, let be a vector in the set (4.1). Then under coordinatewise logarithm, and . Hence is in the image of . Thus is equal to the set (4.1).

It is left to rewrite the inequalities (4.2) as in the statement of the theorem. We havefor all . Here we use Lemma 1 and the definition of the discrete Fourier transformation matrix. If , then , and hence . If , then by Lemma 2. HenceWe replace by and take both sides of the resulting inequality to the power . Finally, making the substitution gives the desired characterization.

For cyclic, Theorem 1 has been independently proven by Baake and Sumner in the context of circulant matrices (Baake and Sumner 2020, Theorem 5.7). Moreover, they show that every embeddable circulant matrix is circulant embeddable (Baake and Sumner 2020, Corollary 5.2).

It follows from Lemma 3 that if a -Markov matrix P is -embeddable, then there exists a unique -rate matrix Q such that . Indeed, since Q and P have both real eigenvalues and the eigenvalues of P are exponentials of eigenvalues of Q, then the eigenvalues of Q are uniquely determined by the eigenvalues of P. Then the -rate matrix Q is the principal logarithm of P.

The inequalities in Theorem 1 imply . Hence the set of -embeddable matrices for a symmetric group-based model is a relatively closed subset of a connected component of the complement of . A relatively closed subset means here a set that can be written as the intersection of a closed subset of and the connected component of the complement of .

In the rest of the current section and in Sect. 5, we will discuss applications of Theorem 1. We will recover known results about -embeddability and as a novel application characterize embeddability for three group-based models of hachimoji DNA.

Example 4

The CFN model is the group-based model associated to the group . The CFN Markov matrices have the formThe discrete Fourier transform matrix isThe eigenvalues of P are and . By Theorem 1, the Markov matrix P is CFN embeddable if and only if or equivalently . This is equivalent to P satisfying , or equivalently . The result that a general stochastic matrix is embeddable if and only if or goes back to Kingman (1962, Proposition 2). Hence P is CFN embeddable if and only if it is embeddable.

Example 5

Recall that the Kimura 3-parameter model is the group-based model associated to group and a K3P Markov matrix P has the form (2.1). The eigenvalues of P areBy Theorem 1, a Markov matrix P is K3P embeddable if and only ifIn the Kimura 2-parameter model and . We get the conditions for the K2P embeddability by setting in (4.4). Hence a K2P Markov matrix is K2P embeddable if and only ifIn the Jukes–Cantor model and . A JC Markov matrix is JC embeddable if and only ifThe K3P embeddability of a K3P Markov matrix with no repeated eigenvalues is equivalent to the embeddability of the matrix. Similarly, the JC embeddability of a JC Markov matrix is equivalent to the embeddability of the matrix. The same is not true for K2P Markov matrices with exactly two coinciding eigenvalues. See Roca-Lacostena and Fernández-Sánchez (2018, Section 3) for similar computations and further discussion on the model embeddability of K3P, K2P, and JC Markov matrices.

Remark 3

By Kingman (1962, Corollary on page 18), the map from rate matrices to transition matrices is locally homeomorphic except possibly when the rate matrix has a pair of eigenvalues differing by a non-zero multiple of . Since for symmetric group-based models rate matrices are real symmetric, then all their eigenvalues are real and hence the map from rate matrices to transition matrices is a homeomorphism. Therefore the boundaries of embeddable transition matrices of symmetric group-based models are images of the boundaries of the rate matrices. For general Markov model, the boundaries of embeddable transition matrices are characterized in Kingman (Kingman 1962, Propositions 5 and 6).

Corollary 1

A -embeddable transition matrix lies on the boundary of the set of -embeddable transition matrices for a symmetric group-based model if and only if it satisfies at least one of the inequalities in Theorem 1 with equality.

Hachimoji DNA

In this section, we suggest three group-based models for a genetic system with eight building blocks recently introduced by Hoshika et al. (2019), and then characterize model embeddability for the proposed group-based models. The genetic system is called hachimoji DNA. It has four synthetic nucleotides denoted S, B, Z, and P in addition to the standard nucleotides adenine (A), cytosine (C), guanine (G) and thymine (T). Detailed descriptions of the four additional nucleotides are given in Hoshika et al. (2019). If in the standard 4-letter DNA, the purines are A and G and the pyrimidines are C and T, then in the hachimoji system, there are additionally purine analogs P and B, and pyrimidine analogs Z and S. The hydrogen bonds occur between the pairs A-T, C-G, S-B and Z-P.

This DNA genetic system with eight building blocks can reliably form matching base pairs and can be read and translated into RNA. It is mutable without damaging crystal structure which is required for molecular evolution. Hachimoji DNA has potential application in bar-coding, retrievable information storage, and self-assembling nanostructures.

The underlying group we suggest for the hachimoji DNA is , since when restricted to the standard 4-letter DNA it gives the group that is the underlying group for the standard DNA models. We identify the nucleotides with the group elements of as follows:The discrete Fourier transformation matrix of the group is

Hachimoji 7-parameter model

The first model we propose is the analogue of the Kimura 3-parameter model and we will call it the hachimoji 7-parameter (H7P) model. In the hachimoji 7-parameter model, each element of the group maps to a distinct label. Thus the labeling function is trivially -compatible. The H7P rate and transition matrices have the formThe eigenvalues of a H7P Markov matrix areBy Theorem 1, such a matrix is H7P embeddable if and only if all eigenvalues are positive and satisfy

Hachimoji 3-parameter model

The second model we suggest specializes to the Kimura 2-parameter model when restricted to the standard 4-letter DNA. We will call it the hachimoji 3-parameter (H3P) model. We recall that in the Kimura 2-parameter model there are three distinct parameters for the rates of mutation: One parameter for a state remaining unchanged, one parameter for transversion from a purine base to a pyrimidine base or vice versa, and one parameter for transition to the other purine or to the other pyrimidine. We say that two bases are of the same type if they are both standard or synthetic bases. In the hachimoji 3-parameter model, there are the following parameters:The H3P rate and transition matrices have the formThe labeling function of this model corresponds to the partitionwhich is -compatible by Table 1.

a: the probability of a state remaining unchanged.

b: the probability of a transversion from a purine base to a pyrimidine base or vice versa.

c: the probability of a transition to another purine or pyrimidine base of the same type (same type transitions).

d: the probability of a transition to another purine or pyrimidine base of different type (different type transitions).

The eigenvalues of a H3P Markov matrix areBy Theorem 1, a H3P Markov matrix P is H3P embeddable if and only if the eigenvalues of P satisfy

Hachimoji 1-parameter model

The third model we suggest is the analogue of the Jukes–Cantor model and we will refer to it as hachimoji 1-parameter (H1P) model. It is the simplest group-based model associated to the group and it is described by only two distinct parameters for the rates of mutation. The two parameters are for a state remaining the same and a state mutating to any other state. The corresponding labeling function is -compatible by Lemma 6. The H1P rate and transition matrices have the formThe eigenvalues of a H1P Markov matrix are and for . By Theorem 1, such a matrix is H1P embeddable if and only if its eigenvalues satisfy

Remark 4

The same conditions as in (5.6) characterize model embeddability for the general Jukes–Cantor model as defined in Lemma 6. This is also a special instance of a more general result (Baake and Sumner 2020, Corollary 4.7) on equal-input embeddability. If the order of is even, then the notion of general embeddability is equivalent to the notion of model embeddability for the general Jukes–Cantor models by Baake and Sumner (2020, Theorem 4.6).

Volume

In this section we compute the relative volumes of model embeddable Markov matrices within some meaningful subsets of Markov matrices by taking advantage of the characterisation of embeddability in terms of eigenvalues. The aim of this section is to describe how large the different sets of matrices are compared to each other and provide intuition of how restrictive is the hypothesis of homogeneous continuous-time models.

We will focus on the hachimoji models and the generalization of the Jukes–Cantor model. We will use the following notation: Biologically, the subspace of diagonally dominant matrices consists of matrices with probability of not mutating at least as large as the probability of mutating. If a diagonally dominant matrix is embeddable, it has an identifiable rate matrix (Cuthbert 1972, 1973), namely a unique Markov generator, which is crucial for proving the consistency of many phylogenetic reconstruction methods, such as those based on maximum likelihood methods (Casanellas et al. 2020c; Chang 1996). What is more, the set of Markov matrices with positive eigenvalues includes the multiplicative closure of the transition matrices in the continuous-time version of the model (Sumner et al. 2012). We have the inclusions and . The volumes of these spaces are given for the Kimura 3-parameter model in Roca-Lacostena and Fernández-Sánchez (2018, Theorem 4.1), for the Kimura 2-parameter model in Casanellas et al. (2020b, Proposition 5.1] and for the Jukes–Cantor model in Roca-Lacostena and Fernández-Sánchez (2018, Section 4).

is the set of all Markov matrices in a model.

is the subset of matrices in with only positive eigenvalues.

is the subset of diagonally dominant matrices in , i.e. matrices in such that in each row the diagonal entry is greater or equal than the sum of all other entries in the same row.

is the subset of model embeddable transition matrices in .

The subsets and can be described using the parameterization in terms of the entries of the Markov matrix or in terms of their eigenvalues. We parameterize the relevant subsets of Markov matrices in terms of the eigenvalues of the Markov matrices and compute the volumes using these parametrizations. If denotes the bijection from the set of entries of a Markov matrix in a particular model to the set of its eigenvalues and the matrix denotes the Jacobian matrix of the map , then the volume of any subset in the parametrization using entries of a Markov matrix will be times the volume in the parameterization using eigenvalues. Since the determinant of this Jacobian is constant for each of the three models we consider, the relative volumes of the set of model embeddable Markov matrices will not depend on the parameterization chosen.

Proposition 1

For the hachimoji 7-parameter model, consider , , and as subsets of parameterized by , the eigenvalues of a H7P Markov matrix. Then: (i) (ii) (iii)

The entries of a H7P Markov matrix (5.2) are determined by a vector (a, b, c, d, e, f, g, h). The entries of this vector can be expressed in terms of the eigenvalues aswhere K is the discrete Fourier transform matrix (5.1). In terms of the entries or the eigenvalues of a H7P Markov matrix, the relevant subsets in this model are given by: is given by one equation and seven inequalities presented in Sect. 5.1. Expressing all conditions defining , , and in terms of the eigenvalues allows us to compute volumes of these sets using Polymake (Gawrilow and Joswig 2000).

We are not able to compute the volume of the subspace of the H7P embeddable Markov matrices exactly. Instead we estimate the volume using the hit-and-miss Monte Carlo integration method (Hammersley 2013) implemented in Mathematica. Table 2 summarizes the volume for various number of sample points. Table 3 gives relative volumes for the relevant sets.

Table 2

The estimated volume of the set of H7P embeddable matrices using the hit-and-miss Monte Carlo integration with n sample points

n	10^4	10^5	10^6	10^7
V(Δ_me)	0.0015	0.00197	0.001946	0.0019678
V(Δ_me∩Δ_dd)	0.0008	0.00084	0.00085	0.0008271

Table 3

The relative volumes for the hachimoji 7-parameter model

	Δ	Δ_+	Δ_dd
V(·)V(Δ)	1	1754096=0.042724609375	1128=0.0078125
V(Δ_me∩·)V(·)	≈0.00239	≈0.056045	≈0.13388

The volumes of and are estimated using Monte Carlo integration with sample points

Proposition 2

For the hachimoji 3-parameter model, consider , , , and as subsets of parameterized by x, y, z, the eigenvalues of a H3P Markov matrix. Then: (i) ; (ii) ; (iii) ; (iv) ; (v) .

The entries of a H3P Markov matrix as in (5.3) can be expressed in terms of the eigenvalues asExpressing all conditions defining , , , and in terms of x, y, z allows us to use the Integrate command in Mathematica to compute the desired volumes. For we used the numerical integration command NIntegrate.

The sets , and for the hachimoji 3-parameter model are depicted in Fig. 1. The relative volumes of relevant sets are given in Table 4.

Fig. 1

The sets , and for the hachimoji 3-parameter model. The sets and are polytopes; the set is a semialgebraic set

Table 4

The relative volumes for the hachimoji 3-parameter model

	Δ	Δ_+	Δ_dd
V(·)V(Δ)	1	2164=0.328125	18=0.125
V(Δ_me∩·)V(·)	14=0.25	1621≈0.76190	≈0.82040

Finally, we discuss the generalization of the Jukes–Cantor model which includes the hachimoji 1-parameter model. Let be a finite abelian group of order n and a labeling function such that and and for . In Lemma 6 we proved that L is a -compatible labeling. In the general Jukes–Cantor model, the transition matrix P corresponding to this labeling is has the formSince P is a Markov matrix, then , and thus P is parameterized by b.

Proposition 3

For the general Jukes–Cantor model, consider , , , as subsets of parameterized by b, the off-diagonal element of the Markov matrix. Then: (i) ; (ii) ; (iii) ; (iv) ; (v) .

The Markov matrix P has eigenvalues 1 with multiplicity 1 and with multiplicity . Hence

(i) .

(ii) .

(iii) .

(iv) By Remark 4, a Markov matrix is general Jukes–Cantor embeddable if and only if the eigenvalue satisfies . Since necessarily holds for any Markov matrix, we have .

(v) Since , then .

The relative volumes of relevant sets for the general Jukes–Cantor model are presented in Table 5. Proposition 3 gives for the hachimoji 1-parameter model (i) ; (ii) ; (iii); (iv) ; (v).

Table 5

The relative volumes for the general Jukes–Cantor model

	Δ	Δ_+	Δ_dd
V(·)V(Δ)	1	n-1n	12
V(Δ_me∩·)V(·)	n-1n	1	1

Conclusion

When modelling sequence evolution we often adopt several simplifying assumptions, which make the statistical problems tractable. The commonly used Markov models depend on the assumption that sites evolve independently following a Markov process. The Markov chain is often assumed to be homogeneous continuous-time, that is the transition probabilities are independent of the time. This means that the instantaneous rates of substitution at any time are fixed and usually displayed as the entries of rate matrices. If an evolutionary process is not homogeneous, then one can multiply transition matrices of short homogeneous processes. The resulting matrix is not necessarily embeddable, but if it is, then the inhomogeneous process can be approximated by a homogeneous one.

In this paper we provide necessary and sufficient conditions for model-embeddability of symmetric group-based substitution models, which include the well known Cavender–Farris–Neyman, Jukes–Cantor, Kimura-2 and Kimura-3 parameter models for DNA. We fully characterize those embeddable stochastic matrices following a symmetric group-based model structure whose Markov generators also satisfy the constraints of the model, which we refer to as model embeddability.

A novel application of our main result is the characterization of model embeddability for three group-based models for the hachimoji DNA, a synthetic genetic system with eight building blocks. For these models we also compute the relevant volumes of model embeddable matrices within other relevant sets of Markov matrices. These computations show how restrictive is the hypothesis of a particular hachimoji time-continuous group-based model.

In this article we have considered symmetric group-based models. The importance of the symmetricity assumption is that it guarantees that the eigenvalues of rate and transition matrices of a group-based model are real. We use this property in the proof of Theorem 1. A future research question is to explore whether this approach can be extended to group-based models that are not symmetric.

Table 6

Friendly labelings for abelian group of order

n	Group	Friendly labelings
2	Z_2	{{0,1}},{{0},{1}}
3	Z_3	{{0,1,2}},{{0},{1,2}},{{0},{1},{2}}
4	Z_4	{{0,1,2,3}},{{0},{1,2,3}},{{0,2},{1,3}},{{0},{1,3},{2}},{{0},{1},{2},{3}}
4	Z_2×Z_2	{{(0,0),(0,1),(1,0),(1,1)}},{{(0,0)},{(0,1),(1,0),(1,1)}}, {{(0,0),(0,1)},{(1,0),(1,1)}},{{(0,0),(1,1)},{(0,1),(1,0)}}, {{(0,0),(1,0)},{(0,1),(1,1)}},{{(0,0)},{(0,1)},{(1,0),(1,1)}}, {{(0,0)},{(0,1),(1,0)},{(1,1)}},{{(0,0)},{(0,1),(1,1)},{(1,0)}}, {{(0,0)},{(0,1)},{(1,0)},{(1,1)}}
5	Z_5	{{0,1,2,3,4}},{{0},{1,2,3,4}},{{0},{1,4},{2,3}},{{0},{1},{2},{3},{4}}
6	Z_2×Z_3	{{(0,0),(0,1),(0,2),(1,0),(1,1),(1,2)}},{{(0,0)},{(0,1),(0,2),(1,0),(1,1),(1,2)}}, {{(0,0),(0,1),(0,2)},{(1,0),(1,1),(1,2)}},{{(0,0),(1,0)},{(0,1),(0,2),(1,1),(1,2)}}, {{(0,0)},{(0,1),(0,2)},{(1,0),(1,1),(1,2)}},{{(0,0)},{(0,1),(0,2),(1,1),(1,2)}, {(1,0)}},{{(0,0),(1,0)},{(0,1),(1,1)},{(0,2),(1,2)}},{{(0,0)},{(0,1)},{(0,2)},{(1,0),(1,1),(1,2)}}, {{(0,0)},{(0,1),(0,2)},{(1,0)},{(1,1),(1,2)}},{{(0,0)},{(0,1),(1,1)},{(0,2),(1,2)},{(1,0)}}, {{(0,0)},{(0,1)},{(0,2)},{(1,0)},{(1,1)},{(1,2)}}
7	Z_7	{{0,1,2,3,4,5,6}},{{0},{1,2,3,4,5,6}},{{0},{1,2,4},{3,5,6}},{{0},{1,6},{2,5},{3,4}},{{0},{1},{2},{3},{4},{5},{6}}
8	Z_8	{{0,1,2,3,4,5,6,7}},{{0},{1,2,3,4,5,6,7},{{0,4},{1,2,3,5,6,7}}, {{0,2,4,6},{1,3,5,7}},{{0},{1,2,3,5,6,7},{4}},{{0},{1,2,6,7},{3,4,5}}, {{0},{1,4,7},{2,3,5,6}},{{0},{1,3,5,7},{2,4,6}},{{0,4},{1,3,5,7},{2,6}}, {{0},{1,3,5,7},{2,6},{4}},{{0,4},{1,5},{2,6},{3,7}},{{0},{1,3,5,7},{2},{4},{6}},{{0},{1,3}{2,6},{4},{5,7}},{{0},{1,7},{2,6},{3,5},{4}},{{0},{1,5},{2,6},{3,7},{4}},{{0},{1,5},{2},{3,7},{4},{6}},{{0},{1},{2},{3},{4},{5},{6},{7}}