Embeddability and rate identifiability of Kimura 2-parameter matrices

TL;DR

This paper fully characterizes embeddable and rate-identifiable Markov matrices within the K80 nucleotide substitution model, revealing sets with non-identifiable rates and implications for phylogenetic parameter estimation.

Contribution

It provides a complete characterization of embeddability and rate identifiability for 4x4 K80 matrices, including open subsets with non-identifiable rates.

Findings

Identified the set of embeddable K80 matrices

Described open subset with non-identifiable rates

Computed relative volumes of embeddable matrices

Abstract

Deciding whether a Markov matrix is embeddable (i.e. can be written as the exponential of a rate matrix) is an open problem even for $4\times 4$ matrices. We study the embedding problem and rate identifiability for the K80 model of nucleotide substitution. For these $4\times 4$ matrices, we fully characterize the set of embeddable K80 Markov matrices and the set of embeddable matrices for which rates are identifiable. In particular, we describe an open subset of embeddable matrices with non-identifiable rates. This set contains matrices with positive eigenvalues and also diagonal largest in column matrices, which might lead to consequences in parameter estimation in phylogenetics. Finally, we compute the relative volumes of embeddable K80 matrices and of embeddable matrices with identifiable rates. This study concludes the embedding problem for the more general model K81 and its submodels, which had been initiated by the last two authors in a separate work.

Tables

Table 1. Table 1 : Relative volumes of the spaces labelling the columns of the table (which are referred to as “ ⋅ ⋅ \,\cdot\, ” in each row) within the set of K80 Markov matrices (first row), relative volumes of embeddable matrices within those spaces (second row) and relative volume of those of them with identifiable rates (third row). The values are rounded to the 10th decimal.

	$\Delta$	${\Delta }_{\mathrm{DLC}}$	${\Delta }_{+}$	${\Delta }_{\mathrm{dd}}$
$\frac{V(\cdot )}{V\left(\Delta \right)}$	1	$\frac{5}{12}\approx 0.4166666667$	0.375	0.25
$\frac{V({\Delta }_{\mathrm{emb}}\cap \cdot )}{V(\cdot )}$	$\frac{1+e^{-3\pi }}{3}\approx$ 0.3477379727	0.8	$\frac{8}{9}\approx 0.8888888889$	$\frac{14-8\sqrt{2}}{3}\approx$ 0.895430501
$\frac{V({\Delta }_{\mathrm{idf}}\cap \cdot )}{V({\Delta }_{\mathrm{emb}}\cap \cdot )}$	$\frac{1-e^{-6\pi }}{1+e^{-3\pi }}\approx$ 0.9999193000	$1-e^{-6\pi }\approx$ 0.9999999935	$1-e^{-6\pi }\approx$ 0.9999999935	1

Equations

$$\left(\begin{array}{cccc}a& b& c& d\\ b& a& d& c\\ c& d& a& b\\ d& c& b& a\end{array}\right),$$

$$M=\left(\begin{array}{cccc}a& b& c& d\\ b& a& d& c\\ c& d& a& b\\ d& c& b& a\end{array}\right).$$

$$S:=\left(\begin{array}{cccc}1& 1& 1& 1\\ 1& 1& -1& -1\\ 1& -1& 1& -1\\ 1& -1& -1& 1\end{array}\right).$$

$$\begin{array}{ccc}a+b+c+d& & x:=a+b-c-d\\ y:=a-b+c-d& & z:=a-b-c+d.\end{array}$$

$$\begin{array}{cccc}\phi :& \text{Entries}& ⟶& \text{Eigenvalues}\\ & (b,c)& ⟼& (1-4c,1-2b-2c)\end{array}\begin{array}{cccc}{\phi }^{-1}:& \text{Eigenvalues}& ⟶& \text{Entries}\\ & (x,y)& ⟼& (\frac{1+x-2y}{4},\frac{1-x}{4})\end{array}$$

$$x,y,z>0,x\ge yz,y\ge xz,z\ge xy.$$

$$Q=S\;U\;\mathrm{diag}\big{(}\log_{k_{1}}(1),\log_{k_{2}}(x),\log_{k_{3}}(y),\log_{k_{4}}(y)\big{)}\;U^{-1}\;S^{-1}$$

$$S^{-1}QS$$

$$S^{-1}QS=\footnotesize\begin{pmatrix}1&0&0&0\\ 0&1&0&0\\ 0&0&\alpha+\beta i&\alpha-\beta i\\ 0&0&\gamma+\delta i&\gamma-\delta i\\ \end{pmatrix}\begin{pmatrix}\;\mathrm{diag}\big{(}0,\log(x)\big{)}\;&0\\ 0&\mathrm{diag}\big{(}\log_{k_{3}}(y),\overline{\log_{k_{3}}(y)}\big{)}\\ \end{pmatrix}\footnotesize\begin{pmatrix}1&0&0&0\\ 0&1&0&0\\ 0&0&\alpha+\beta i&\alpha-\beta i\\ 0&0&\gamma+\delta i&\gamma-\delta i\\ \end{pmatrix}^{-1}.$$

$$R=\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\\ 0& 0& 1& 1\\ 0& 0& i& -i\end{array}\right),$$

$$L_{0}=S\;U\;\mathrm{diag}\big{(}0,\log(x),\log|y|,\log|y|\big{)}\;U^{-1}\;S^{-1}\,.$$

$$L_1(A)=S(U{R}^{-1})\mathrm{diag}\left(\left(\begin{array}{cc}0& 0\\ 0& 0\end{array}\right),\left(\begin{array}{cc}0& 1\\ -1& 0\end{array}\right)\right)(U{R}^{-1}{)}^{-1}{S}^{-1}=(SU)\mathrm{diag}(0,0,i,-i)(U^{-1}{S}^{-1}).$$

$$Q(k,A)=(SU)\mathrm{diag}(\log (1),\log (x),{\log }_k(y),{\log }_l(y))(U^{-1}{S}^{-1})$$

$$Q(k,\alpha ,\beta ):=Q(k,A)\text{ with }k\in \mathbb{Z}\text{ and }A=\left(\begin{array}{cc}1& 0\\ \alpha & \beta \end{array}\right)\text{ for some }\alpha \in \mathbb{R},\beta \in {\mathbb{R}}_{>0}.$$

$$Q(k,A)=\{\begin{array}{ll}Q\left(k,\frac{a_{11}{a}_{21}+a_{12}{a}_{22}}{a_{11}^2+a_{12}^2},\frac{\det (A)}{a_{11}^2+a_{12}^2}\right)& \text{ if }\det (A)>0,\text{vspace}\ast 5mm\\ Q\left(-k-\frac{{\theta }_y}{\pi },\frac{a_{11}{a}_{21}+a_{12}{a}_{22}}{a_{11}^2+a_{12}^2},-\frac{\det (A)}{a_{11}^2+a_{12}^2}\right)& \text{ if }\det (A)<0.\end{array}$$

$$S\mathrm{diag}(I{d}_2,A)\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\\ 0& 0& 1& 1\\ 0& 0& i& -i\end{array}\right)=\left(\begin{array}{cccc}1& 1& 1& 1\\ 1& 1& -1& -1\\ 1& -1& 1& -1\\ 1& -1& -1& 1\end{array}\right)\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\\ 0& 0& 1& 0\\ 0& 0& \alpha & \beta \end{array}\right)\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\\ 0& 0& 1& 1\\ 0& 0& i& -i\end{array}\right).$$

$$\displaystyle\log|y|\pm\big{(}2\pi k+\theta_{y}\big{)}i,$$

$$u_{+}(\alpha ,\beta )$$

$$u_{-}(\alpha ,\beta )$$

$$\displaystyle\log|y|\pm\big{(}2\pi k+\theta_{y}\big{)}i=\log|y|\pm\big{(}2\pi k^{\prime}+\theta_{y}\big{)}i.$$

$$\log (x)-2\log |y|$$

$$-\log (x)$$

$$L_0=1/4\left(\begin{array}{cccc}-\lambda -2\mu & -\lambda +2\mu & \lambda & \lambda \\ -\lambda +2\mu & -\lambda -2\mu & \lambda & \lambda \\ \lambda & \lambda & -\lambda -2\mu & -\lambda +2\mu \\ \lambda & \lambda & -\lambda +2\mu & -\lambda -2\mu \end{array}\right).$$

$$L_{1}(A)=\frac{1}{4\beta}\small\begin{pmatrix}1-\alpha^{2}-\beta^{2}&-(1-\alpha^{2}-\beta^{2})&-\big{(}(1+\alpha)^{2}+\beta^{2}\big{)}&(1+\alpha)^{2}+\beta^{2}\vspace*{3mm}\\ -(1-\alpha^{2}-\beta^{2})&1-\alpha^{2}-\beta^{2}&(1+\alpha)^{2}+\beta^{2}&-\big{(}(1+\alpha)^{2}+\beta^{2}\big{)}\vspace*{3mm}\\ (1-\alpha)^{2}+\beta^{2}&-\big{(}(1-\alpha)^{2}+\beta^{2}\big{)}&-(1-\alpha^{2}-\beta^{2})&1-\alpha^{2}-\beta^{2}\vspace*{3mm}\\ -\big{(}(1-\alpha)^{2}+\beta^{2}\big{)}&(1-\alpha)^{2}+\beta^{2}&1-\alpha^{2}-\beta^{2}&-(1-\alpha^{2}-\beta^{2})\\ \end{pmatrix}\,.$$

$$\displaystyle-\lambda+2\mu\pm\big{(}2\pi k+\theta_{y}\big{)}\frac{1-\alpha^{2}-\beta^{2}}{\beta}$$

$$\displaystyle\lambda\pm\big{(}2\pi k+\theta_{y}\big{)}\frac{(1+\alpha)^{2}+\beta^{2}}{\beta}$$

$$\displaystyle\lambda\pm\big{(}2\pi k+\theta_{y}\big{)}\frac{(1-\alpha)^{2}+\beta^{2}}{\beta}$$

$$I_k=\{\begin{array}{ll}⟨-k,k⟩& \text{if }y>0\\ ⟨-k-1,k⟩& \text{if }y<0\end{array}$$

$$\log (x)-2\log (y)$$

$$-\log (x)$$

$$Q(0,0,1)=S\left(\begin{array}{cccc}0& 0& 0& 0\\ 0& \log (x)& 0& 0\\ 0& 0& \log |y|& -{\theta }_y\\ 0& 0& {\theta }_y& \log |y|\end{array}\right)S^{-1}.$$

$$M=\left(\begin{array}{cccc}0.2509345932& 0.2490671504& 0.2499991282& 0.2499991282\\ 0.2490671504& 0.2509345932& 0.2499991282& 0.2499991282\\ 0.2499991282& 0.2499991282& 0.2509345932& 0.2490671504\\ 0.2499991282& 0.2499991282& 0.2490671504& 0.2509345932\end{array}\right).$$

Full text

Keywords. Nucleotide substitution model; Markov matrix; Markov generator; matrix logarithm; embedding problem; rate identifiability.

1 Introduction

Modeling molecular substitutions is the first and arguably the most fundamental step in phylogenetics and it is where the crucial hypotheses that should allow us to reconstruct the evolutionary history are to be assumed. There are several approaches one might take and it is natural to ask which conditions must be required for evolutionary models to fit the real evolutionary processes. For instance, under a Markov process, which are the plausible structures for substitution matrices in a given nucleotide or amino acid substitution model? Should any Markov matrix with such structure be considered as a biologically realistic substitution matrix? Concerned with these and similar questions, we address the study of the embedding problem for the matrices of nucleotide substitution models. The final aim of this problem is to characterize those Markov matrices that are consistent with a homogeneous continuous-time approach of evolution. In other words, the aim is to decide whether a given Markov matrix in the model can be written as the exponential of some rate matrix. In this case, the Markov matrix is said to be embeddable (or to have a continuous realization as in the book by Steel, 2016), and the rate matrix is called a Markov generator of the Markov matrix.

The motivation and applications of the embedding problem are diverse and include economics, social sciences, and molecular evolution (see Israel et al.,, 2001; Singer and Spilerman,, 1976; Jia et al.,, 2014; Jia,, 2016; Verbyla et al.,, 2013, for example). The problem has been addressed in a number of papers (e.g. Culver,, 1966; Cuthbert,, 1972, 1973; Davies,, 2010; Goodman,, 1970; Guerry,, 2013; Van-Brunt,, 2018; Guerry,, 2019), and it has been solved for $2\times 2$ and $3\times 3$ matrices by Israel et al., (2001); Goodman, (1970); Guerry, (2019). However, it is far from being solved for larger matrices with full generality. Because our motivation and interest arises from the study of mathematical properties of nucleotide substitution models, we focus on the particular case of $4\times 4$ matrices and, more specifically, on the Kimura models of nucleotide substitution K80 and K81. These models were proposed by Kimura, (1980, 1981) in the celebrated papers published almost 40 years ago. The main motivation was to include different parameters for different nucleotide substitution types, which according to biological data are not equally likely. The substitution matrices of nucleotide replacements for the K81 model are Markov matrices with the following structure

[TABLE]

where the rows and columns are labelled by nucleotides adenine, guanine, cytosine and thymine (in this order) and the entry $(i,j)$ is the conditional probability that nucleotide $i$ is replaced by nucleotide $j$. The K80 model is the submodel obtained when imposing that $c=d$. In this case the parameter $c=d$ corresponds to transversions, which are substitutions from purines (adenine and guanine) to pyrimidines (cytosine and thymine) and viceversa, and parameter $b$ corresponds to transitions (substitutions within purines or within pyrimidines). When $b=c=d$, the resulting submodel is the Jukes-Cantor model JC69 introduced by Jukes and Cantor, (1969).

In the phylogenetic setting, the embedding problem has been recently studied by Jia, (2016) and Roca-Lacostena and Fernández-Sánchez, (2018). In the first paper, the author addresses the embedding problem restricted to time-reversible Markov generators. In the second, the authors characterize the embeddability of K81 matrices with different eigenvalues. In a parallel work, Kosta and Kubjas, (2017) characterize those embeddable matrices in a group-based model whose Markov generators satisfy the constraints of the model; this is usually referred to as model embeddability. Although the problem of model embeddability is natural from a continuous-time approach, in this paper we do not impose time-reversibility or any model restriction on rates, as we are more interested in determining whether a Markov process has a homogeneous continuous-time realization or not (independently of the structure of the instantaneous rate matrix). On the other hand, Kimura, (1981) initiates his argument in terms of differences of nucleotides between two homologous sequences, so his model for transition and transversion parameters is firstly justified in terms of probabilities rather than instantaneous mutation rates. As shown here, rates do not necessarily satisfy the same symmetries as probabilities do (see Example 4.3).

Embeddability has not been yet fully characterized for K81 matrices with a negative eigenvalue of even multiplicity and, in particular, it has not been solved for matrices within the K80 model. In the present work we address the study of the repeated eigenvalues case to fill this gap and give a complete answer for the K80 model. Moreover, together with the results by Roca-Lacostena and Fernández-Sánchez, (2018), our results fully solve the embedding problem for the K81 model. More precisely, the results on the embeddability of K81 matrices can be summarized as follows (see also Remark 3.12):

• (1)

for generic K81 Markov matrices (where $b,c,d$ are different), embeddability holds if and only if the principal logarithm $Log(M)$ is a rate matrix;

• (2)

for K81 Markov matrices with exaclty two equal off-diagonal entries, embeddability is characterized by the theorem below (and its analogous versions under permutations of rows and columns, i.e., if $b=d$ one has to exchange the roles of $b$ and $c$, and if $b=c$, then $b$ should be replaced by $d$.);

• (3)

for the JC69 model, embeddability holds if and only if the determinant is positive.

Moreover, in this paper, we go further in the study of the K80 model and derive a criterion to determine whether the mutation rates of embeddable matrices are identifiable and, if not, we determine whether there is a finite number of Markov generators. One says that the rates of a Markov matrix $M$ are identifiable there is a unique Markov generator for $M$. The concept of non-identifiability of rates was called nonunique mapping in the paper by Kaehler et al., (2015).

All together, the main results of this work for the K80 model are summarized in the following result:

Theorem 1.1**.**

For any K80 Markov matrix $M=\footnotesize\begin{pmatrix}a&b&c&c\\ b&a&c&c\\ c&c&a&b\\ c&c&b&a\end{pmatrix}$ with $b\neq c$, the following holds:

1. (a)

If ${2}c=1-2b$, then $M$ is not embeddable. 2. (b)

If ${2}c<1-2b$, $M$ is embeddable if and only if $c\leq\sqrt{b}-b$. In this case, $Log(M)$ is a Markov generator and

1. i)

if $c<\frac{1}{4}-\frac{e^{-4\pi}}{4}$ then the rates of $M$ are identifiable, 2. ii)

if $c=\frac{1}{4}-\frac{e^{-4\pi}}{4}$ then $M$ has exactly $3$ Markov generators, 3. iii)

if $c>\frac{1}{4}-\frac{e^{-4\pi}}{4}$ then $M$ has infinitely many Markov generators. 3. (c)

If ${2}c>1-2b$, $M$ is embeddable if and only if $\frac{1}{4}-\frac{e^{-2\pi}}{4}\leq c\leq\sqrt{b}-b$. In this case, the rates of $M$ are not identifiable and

1. i)

if $c=\frac{1}{4}-\frac{e^{-2\pi}}{4}$ then $M$ has exactly $2$ Markov generators. 2. ii)

if $c>\frac{1}{4}-\frac{e^{-2\pi}}{4}$ then $M$ has infinitely many Markov generators.

Case (b) above corresponds to the case of positive eigenvalues, for which we have that $M$ is embeddable if and only if its principal logarithm is a rate matrix (Corollary 3.8). This result is in accordance with the well known case of different and real eigenvalues (see Culver, 1966 and (1) above). Besides this, the theorem also provides a full description of those embeddable K80 matrices with positive eigenvalues and non-identifiable rates (see also Roca-Lacostena and Fernández-Sánchez,, 2018, Example 4.2 and Section 3). Among these K80 matrices, there are some diagonal largest in column matrices (DLC), which are of relevance in parameter estimation in phylogenetics (see Remark 4.9 and the Discussion). As far as we are aware, this is the first result in this direction.

From the results obtained here, we recover the fact that among all possible Markov generators of an embeddable K80 matrix at most one keeps the K80 structure, which in turn coincides with the principal logarithm of the matrix (see Roca-Lacostena and Fernández-Sánchez,, 2018; Jia,, 2016).

Another goal of the paper is to quantify the difference between restricting to the homogeneous continuous-time models or considering any Markov matrix within the K80 model. In this direction, we compute the relative volume of embeddable K80 matrices inside the whole set of K80 Markov matrices. It turns out that embeddable matrices only account for about 35% of all the substitution matrices. Similar computations exhibit that, although embeddable matrices with non-identifiable rates describe a set of positive measure within the K80 model (containing the embeddable matrices in case $(c)$ above), 99.99% of embeddable K80 matrices have identifiable rates, so the non-identifiability of rates is not representative of the general situation (see Table 1 for the precise figures).

The organization of the paper is as follows. In Section 2 we introduce K80 matrices, state the embedding and rate-identifiability problem with precision and recall some known results that are relevant for posterior work. In Section 3 we characterize embeddability for K80 matrices both in terms of their eigenvalues (Corollary 3.9) and their entries (Corollary 3.11), which together with the results already known for K81 matrices do fully solve the embedding problem for the K81 model. Later in Section 4 we solve the rate identifiability problem by providing sufficient and necessary conditions for the rates to be identifiable (Proposition 4.1) and characterize rate identifiability in terms of the eigenvalues or the entries of the Markov matrix (Theorem 4.5 and Corollary 4.8). In Section 5 we compute the relative volumes of embeddable and rate identifiable matrices within the K80 model and some biologically relevant subsets (see Table 1). Finally, in Section 6 we discuss the implications and connections with other papers and possibilities for future work.

2 Preliminaries

In this section we introduce the embedding problem for Markov matrices, the nucleotide substitution model we work with, and a few known results related to the embedding problem needed for the sequel.

Definition 2.1**.**

A matrix $M\in M_{n}(\mathbb{R})$ is said to be a Markov matrix (or substitution matrix) if its entries are positive or zero and its rows sum up to one. Similarly, $Q\in M_{n}(\mathbb{R})$ is said to be a rate matrix if its rows sum to zero and its off-diagonal entries are non-negative.

If $Q$ is a rate matrix, it is well-known that $e^{tQ}=\sum_{n\geq 0}\frac{t^{n}Q^{n}}{n!}$ is a Markov matrix for all $t\geq 0$, thus $Q$ is said to be a Markov generator for $M=e^{Q}$ (Davies,, 2010). However, not every Markov matrix can be obtained this way. A Markov matrix $M$ is said to be embeddable if $M=e^{Q}$ for some rate matrix $Q$. Characterizing which Markov matrices are embeddable is known as the embedding problem. In terms of Markov processes, a Markov matrix is embeddable if it has a realization as a homogeneous continuous-time Markov chain ($M=e^{tQ}$ for some $t>0$ and some rate matrix $Q$).

We say that $Q\in M_{n}(\mathbb{C})$ is a logarithm of a matrix $M$ if $e^{Q}=M$. As in the case of complex numbers, matrix logarithms are not unique (Gantmacher,, 1959). Hence an embeddable matrix may admit more than one Markov generator. If an embeddable matrix has only one Markov generator we say that its rates are identifiable. The identifiability problem consists on deciding whether the rates are identifiable or not.

Definition 2.2**.**

Given $z\in\mathbb{C}\setminus\{0\}$, let $\theta_{z}\in(-\pi,\pi]$ be its principal argument. For any given $k\in\mathbb{Z}$ we define the k-th determination of its logarithm as $\log_{k}(z)=\log(|z|)+(\theta_{z}+2\pi k)i$. We denote by $\log(z)$ its principal logarithm $\log_{0}(z)$. The principal logarithm of an invertible matrix $M$, denoted as $\mathrm{Log}(M)$, is defined as the unique logarithm such that the imaginary part of its eigenvalues lies in the interval $(-\pi,\pi]$ (Higham,, 2008). If $M$ diagonalizes, $M=P\;\mathrm{diag}(\lambda_{1},\dots,\lambda_{n})\;P^{-1}$, then its principal logarithm can be computed as $\mathrm{Log}(M)=P\;\mathrm{diag}(\log(\lambda_{1}),\dots,\log(\lambda_{n}))\;P^{-1}$.

Our goal in this work is to deal with both the embedding problem and the identifiability problem for Kimura 2-substitution types model (K80 for short) (Kimura,, 1980) which is a submodel of Kimura 3-substitution types model K81 (Kimura,, 1981) and contains the most simple Jukes-Cantor model JC69 (Jukes and Cantor,, 1969).

Definition 2.3**.**

A matrix $M\in M_{4}(\mathbb{R})$ is a K81 matrix if it is of the form

[TABLE]

For ease of reading, we will use the notation $M=K(a,b,c,d)$. A K80 matrix is a K81 matrix with $c=d$. Similarly, a JC69 matrix is a K81 matrix with $b=c=d$.

When a K81 matrix is also a Markov matrix, we speak about K81 Markov matrices (respectively K80 Markov matrices and JC69 Markov matrices). When the rows of a K81 matrix sum to zero and the off-diagonal entries are non-negative, then it is a rate K81 matrix (and analogously for K80 and JC69).

As a variation of the embedding problem one may consider the problem of model embeddability mentioned in the Introduction (which attempts to determine those Markov matrices in a certain model that have a Markov generator preserving the same identities that characterize the model). For instance, one may want to characterize those Markov matrices for which there is a K81 Markov generator. For this model (and its submodels), model embeddability has been characterized in terms of inequalities in the eigenvalues of the Markov matrices (Roca-Lacostena and Fernández-Sánchez,, 2018). The model embeddability for general group based models has been studied by Kosta and Kubjas, (2017). As we explained in the Introduction, we shall not focus in this problem in the present paper.

It is well known that all K81 matrices diagonalize under the following Hadamard matrix (Evans and Speed,, 1993; Hendy and Penny,, 1993):

[TABLE]

Furthermore, it follows from a straightforward computation that a K81 matrix has the following eigenvalues:

[TABLE]

In particular, any K81 matrix has real eigenvalues. If we deal with Markov matrices, then the eigenvalues become $1=a+b+c+d$, $x=1-2c-2d$, $y=1-2b-2d$, and $z=1-2b-2c$. For a K80 Markov matrix the eigenvalues become $1$, $x=1-4c$, and $y=z=1-2b-2c$ . Similarly, the eigenvalues of a JC69 Markov matrix are $1$ and $x=y=z=1-4b$. Note that this allows us to parametrize all Markov K81 matrices by the entries $b,c,d$ ($a=1-b-c-d$) or by their eigenvalues $x,y,z$. In particular, for the K80 model, the bijection between both spaces of parameters is given by:

[TABLE]

Remark 2.4**.**

It is worth noting that if $M$ is a K80 Markov matrix with eigenvalues $1,x,y,y$ then

1. a)

$|x|\leq 1$, $|y|\leq 1$ (by Perron-Frobenius theorem). 2. b)

If $y\neq x$, then $y\neq 1$. Indeed, the equality (2) shows that the eigenvalues of $M=K(1-b-2c,b,c,c)$ are $x=1-4{c}$ and $y=1-2b-2c$. Since $b,c\geq 0$, if $y=1$ it follows that $b=c=0$ and $x=1$, which contradicts $y\neq x$. Therefore we have $y\neq 1$. 3. c)

If $y\neq x$ and $x<0$, then $M$ does not have any real logarithm (see Culver,, 1966). 4. d)

If $x$ or $y$ are zero, then $M$ is not embeddable because it has zero determinant (note that $\det e^{Q}=e^{\mathrm{tr}(Q)}>0$). 5. e)

If $x=y$ then $M$ is a JC69 matrix. In this case, $M$ is embeddable if and only if $x>0$ (Roca-Lacostena and Fernández-Sánchez,, 2018).

It is known that if a matrix has determinant close to $1$ then the principal logarithm is the only possible real logarithm of that matrix (Israel et al.,, 2001; Cuthbert,, 1972, 1973; Singer and Spilerman,, 1976). The same holds for matrices with distinct real eigenvalues (Culver,, 1966). The next theorem provides sufficient and necessary conditions for the principal logarithm of a K81 Markov matrix to be a rate matrix (Roca-Lacostena and Fernández-Sánchez,, 2018, see also example 4.5 by Kosta and Kubjas, 2017). In particular, the result below solves the embedding and identifiability problems in the K81 model except for those matrices with repeated eigenvalues. At the same time, this result solves the model embeddability problem for the K81 model and its submodels.

Theorem 2.5** (Roca-Lacostena and Fernández-Sánchez, 2018, Corollary 3.5).**

Let $M$ be a K81 Markov matrix with eigenvalues $1,x,y,z$. Then,

1. i)

$\mathrm{Log}(M)=S\;\mathrm{diag}(0,\log(x),\log(y),\log(z))\;S^{-1}$* is a K81 matrix. Furthermore, it is the only logarithm of $M$ that is itself a K81 matrix.* 2. ii)

$\mathrm{Log}(M)$* is a rate matrix if and only if*

[TABLE]

If a Markov matrix has a repeated eigenvalue (e.g. K80 and JC69 matrices), then it may have infinitely many real logarithms (Culver,, 1966) and hence we must check if any of them is a rate matrix before deciding that such a matrix is not embeddable. Nonetheless, it is known that any logarithm (including non-real logarithms) of a given matrix $M$ can be obtained as $\mathrm{Log}(M)+L$ where $L$ is one out of infinite logarithms of the identity that commute with $M$ (Higham,, 2008).

The following notation will be used throughout the paper.

Note. $Id_{n}$ denotes the identity matrix of order $n$. We write $GL_{n}(\mathbb{K})$ for the space of $n\times n$ invertible matrices with entries in $\mathbb{K}=\mathbb{R}$ or $\mathbb{C}$. Given a matrix $M$, we define the commutant of $M$, $Comm^{*}(M)$, as the set of invertible complex matrices that commute with $M$.

Remark 2.6**.**

If $D$ is a diagonal matrix, $D=\mathrm{diag}(\overbrace{\lambda_{1},\dots,\lambda_{1}}^{m_{1}},\overbrace{\lambda_{2},\dots,\lambda_{2}}^{m_{2}},\dots,\overbrace{\lambda_{n},\dots,\lambda_{n}}^{m_{n}})$, then $Comm^{*}(D)$ consists on all the block-diagonal matrices whose blocks are taken from the corresponding $GL_{m_{i}}(\mathbb{C})$. In particular, if all $\lambda_{i}$ are different, $Comm^{*}(D)$ is the set of invertible diagonal matrices. Moreover, if $M=PM^{\prime}P^{-1}$, then $Comm^{*}(M)$ is formed by all matrices that arise as $P\,U\,P^{-1}$ with $U\in Comm^{*}(M^{\prime})$.

Using that all K80 matrices can be diagonalized by $S$, immediate application of Theorem 1.27 by Higham, (2008) allows us to enumerate all the logarithms (real or not) of any K80 Markov matrix:

Theorem 2.7**.**

Given a K80 Markov matrix $M$ with eigenvalues $1$, $x$, $y$, $y$ all (complex) solutions to $exp(Q)=M$ are given by:

[TABLE]

where $k_{i}\in\mathbb{Z}$, $S$ is defined in (1) and U\in Comm^{*}\big{(}\mathrm{diag}(1,x,y,y)\big{)}.

The previous theorem provides the eigendecomposition of any complex logarithm of a K80 Markov $M$ while making explicit the connection with the commutant of $M$ (see Remark 2.6). In the next section, we will consider only real logarithms, so we will need to pay special attention to the corresponding restrictions on their eigenvalues and eigenvectors.

3 Embeddability of K80 Markov matrices

In this section we characterize all the real logarithms with rows summing to zero of K80 Markov matrices and, as a consequence, we are able to provide sufficient and necessary conditions for such a matrix to be embeddable. These conditions are given in terms of some inequalities involving either the eigenvalues of the matrix (Corollary 3.9) or its entries (Corollary 3.11).

Theorem 2.7 allows us to compute all the logarithms of any given K80 matrix by using the principal logarithm and logarithms of the identity. The main issue with the description given by this theorem is that we are interested only in logarithms that satisfy rate matrices constraints. The forthcoming Proposition 3.2 characterizes all the real logarithms of a K81 matrix whose rows sum to 0. To this end we introduce the following matrices:

Definition 3.1**.**

Let $M$ be a K80 Markov matrix with eigenvalues $1$, $x$, $y$, $y$ satisfying $x>0$, $y\neq x$ and $y\neq 0$. Given $k\in\mathbb{Z}$ and $A\in GL_{2}(\mathbb{R})$ we introduce the following notation:

• •

L_{0}:=S\;\mathrm{diag}\big{(}0,\log(x),\log|y|,\log|y|\big{)}\;S^{-1}.

• •

L_{1}(A):=\big{(}S\;\mathrm{diag}(Id_{2},A)\big{)}\;\mathrm{diag}\left(\footnotesize\begin{pmatrix}0&0\\ 0&0\\ \end{pmatrix},\footnotesize\begin{pmatrix}0&1\\ -1&0\\ \end{pmatrix}\right)\big{(}\mathrm{diag}\big{(}Id_{2},A^{-1}\big{)}\;S^{-1}\big{)}.

• •

Q(k,A):=L_{0}+\big{(}2\pi k+\theta_{y}\big{)}L_{1}(A).

Although these matrices depend on $M$, for ease of reading we decided not to reflect it in the notation. Note that $\theta_{y}$ is either [math] or $\pi$ because $y$ is a real number (see (2)). Also note that if $y>0$, then $L_{0}=Log(M)$ and $Q(0,A)=Log(M)$ for all $A\in GL_{2}(\mathbb{R})$.

Proposition 3.2**.**

Let $M$ be an invertible K80 Markov matrix with eigenvalues $1$, $x$, $y$, $y$ satisfying $x>0$ and $y\neq x$. Then for any matrix $Q\in M_{4}(\mathbb{R})$ the following are equivalent:

1. i)

$Q$* is a real logarithm of $M$ whose rows sum to 0.* 2. ii)

$Q=Q(k,A)$* for some $k\in\mathbb{Z}$ and $A\in GL_{2}(\mathbb{R})$.*

It is worth pointing out that the previous definition and proposition can be generalized to $4\times 4$ diagonalizable matrices. We do not state them in this full generality in order to make the notation and the proof more readable. The proof of this result relies on Theorem 2.7 and the study of the commutant of K80 matrices.

Proof.

Note that Remark 2.4 shows that $y\neq 1$.

$i)\Rightarrow ii)$ Assume that $Q$ is a real logarithm of $M$ whose rows sum to 0. According to Theorem 2.7 it holds that

[TABLE]

for some $k_{i}\in\mathbb{Z}$ and some $U\in GL_{4}(\mathbb{C})$ that commutes with \mathrm{diag}\big{(}1,x,y,y\big{)}.

Since the rows of $Q$ sum to [math], we have that $(1,1,1,1)$ is an eigenvector of $Q$ with eigenvalue [math]. Moreover, since the matrix $Q$ is real, complex non-real eigenvalues of $Q$ come in conjugate pairs. Therefore, as $y\neq x$ and $y\neq 1$, $\log_{k_{1}}(1)$, $\log_{k_{2}}(x)$ are real and at least one of them is equal to [math]. In particular, $k_{1}=k_{2}=0$.

Noting that $y\neq 1$, $y\neq x$ and using Remark 2.6 we get:

$U=\begin{pmatrix}U_{1}&0\\ 0&U_{2}\\ \end{pmatrix}$ with U_{1}\in Comm^{*}\big{(}\mathrm{diag}(1,x)\big{)} and $U_{2}\in GL_{2}(\mathbb{C})$.

Now, as $U_{1}$ commutes with $\mathrm{diag}(1,x)$, it does also commute with \mathrm{diag}\big{(}\log_{0}(1),\log_{0}(x)\big{)}=\mathrm{diag}\big{(}0,\log(x)\big{)} and hence:

[TABLE]

Since non-real eigenvalues and eigenvectors of $Q$ appear in complex conjugate pairs, either $\log_{k_{3}}(y)$, $\log_{k_{4}}(y)\in\mathbb{R}$ or $\log_{k_{4}}(y)=\overline{\log_{k_{3}}(y)}$ and the vector columns of $U_{2}$, namely $v$ and $w$, must be complex vectors satisfying $w=\lambda\overline{v}$ for some $\lambda\in\mathbb{C}\setminus\{0\}$. In the first case we have that $k_{3}=k_{4}=0$ and $y>0$ so that $U_{2}$ commutes with \mathrm{diag}\big{(}\log_{0}(y),\log_{0}(y)\big{)}. Note that in this case we have $Q=Log(M)$ and coincides with $Q(0,A)$ for any $A$. For the second case we get:

U_{2}=\big{(}v\;w\big{)}=\begin{pmatrix}\alpha+\beta i&\alpha-\beta i\\ \gamma+\delta i&\gamma-\delta i\\ \end{pmatrix}\;\mathrm{diag}(1,\lambda)

for some $\lambda\in\mathbb{C}\setminus\{0\}$, $\alpha,\beta,\gamma,\delta\in\mathbb{R}$ such that $\det(U_{2})={2}\lambda(\beta\gamma-\alpha\delta)i\neq 0$.

Since the matrix $\mathrm{diag}(1,\lambda)$ commutes with any diagonal matrix (and in particular with $\mathrm{diag}(\log_{k_{3}}(y),\overline{\log_{k_{3}}(y)})$, we obtain that (7) is equivalent to the following equality:

[TABLE]

Furthermore, by considering the following invertible matrix

[TABLE]

we have

[TABLE]

Using the distributive property and the fact that the invertible matrix $A:=\begin{pmatrix}\alpha&\beta\\ \gamma&\delta\\ \end{pmatrix}\in GL_{2}(\mathbb{R})$ commutes with $\mathrm{diag}(\log|y|,\log|y|)$ we can rewrite $Q$ as:

S\;\begin{pmatrix}0&0&0&0\\ 0&\log(x)&0&0\\ 0&0&\log|y|&0\\ 0&0&0&\log|y|\\ \end{pmatrix}\;S^{-1}+\big{(}2\pi k_{3}+\theta_{y}\big{)}\;S\;\begin{pmatrix}1&0&0&0\\ 0&1&0&0\\ 0&0&\alpha&\beta\\ 0&0&\gamma&\delta\\ \end{pmatrix}\;\begin{pmatrix}0&0&0&0\\ 0&0&0&0\\ 0&0&0&1\\ 0&0&-1&0\\ \end{pmatrix}\;\begin{pmatrix}1&0&0&0\\ 0&1&0&0\\ 0&0&\alpha&\beta\\ 0&0&\gamma&\delta\\ \end{pmatrix}^{-1}\;S^{-1}.

That is, Q=L_{0}+\big{(}2\pi k_{3}+\theta_{y}\big{)}L_{1}(A) and this concludes this part of the proof.

$ii)\Rightarrow i)$ By definition $Q(k,A)$ is a real matrix. Furthermore, $(1,1,1,1)$ is an eigenvector with eigenvalue [math] for both $L_{0}$ and $L_{1}$ thus the rows of $Q(k,A)$ sum to [math]. A straightforward computation shows that $exp(Q(k,A))=M$ for any $A\in GL_{2}(\mathbb{R})$ and $k\in\mathbb{Z}$. Indeed, note that the matrix $R$ defined in (8) and the matrix $\mathrm{diag}(Id_{2},A)$ both commute with $\mathrm{diag}(a,b,c,c)$ for any $a,b,c\in\mathbb{R}$ so, the matrix defined as $U:=\mathrm{diag}(Id_{2},A)\;R$ commutes with $\mathrm{diag}(0,\log(x),\log|y|,\log|y|)$. In particular we can write $L_{0}$ as

[TABLE]

On the other hand, we have

[TABLE]

Hence,

[TABLE]

with $l=-k$ if $y>0$ and $l=-k-1$ if $y<0$. As $U$ also commutes with $\mathrm{diag}(1,x,y,y)$, the matrix $Q(k,A)$ is one of the matrices listed in Theorem 2.7 and so, it is a logarithm of $M$.

∎

By using the proposition above we obtain a parametrization of all real logarithms with rows summing to 0 of any K80 Markov matrix. Nonetheless, this parametrization is not injective because distinct choices of $A$ and $k$ may produce the same logarithm of $M$. Theorem 3.5 below provides an injective parametrization for these logarithms (other than $\mathrm{Log}(M)$) by considering only matrices of the following form:

[TABLE]

Remark 3.3**.**

We claim that for any $A=(a_{ij})\in GL_{2}(\mathbb{R})$ and any $k\in\mathbb{Z}$, there exist $\alpha,\beta,\tilde{k}$ such that $Q(k,A)=Q(\tilde{k},\alpha,\beta)$. Indeed, one can easily check that $L_{1}(A)$ can be realized as

$L_{1}(A)=L_{1}\left(\begin{pmatrix}1&0\\ \frac{a_{11}a_{21}+a_{12}a_{22}}{a_{11}^{2}+a_{12}^{2}}&\frac{\det(A)}{a_{11}^{2}+a_{12}^{2}}\\ \end{pmatrix}\right)$ and $L_{1}(A)=-L_{1}\left(\begin{pmatrix}1&0\\ \frac{a_{11}a_{21}+a_{12}a_{22}}{a_{11}^{2}+a_{12}^{2}}&-\frac{\det(A)}{a_{11}^{2}+a_{12}^{2}}\\ \end{pmatrix}\right)$.

Since we want the parameter $\beta$ to be positive, we will consider one expression or the other depending on the sign of $\det(A)$. Since $L_{0}$ does not depend on $A$, we get that:

[TABLE]

Note that $\frac{\theta_{y}}{\pi}$ is either 0 (if $y>0$) or 1 (if $y<0$). Summing up, by Proposition 3.2, any real logarithm $Q$ of $M$ with rows summing to [math] can be expressed as $Q=Q(k,\alpha,\beta)$ for some $\alpha,\beta$.

Remark 3.4**.**

Note that the equality (9) obtained in the proof of Proposition 3.2 allows the computation of eigenvalues and eigenvectors for $Q(k,\alpha,\beta)$. Indeed, (9) gives us that the eigenvalues of $Q(k,A)$ are $\log(1),\log(x),\log_{k}(y),\overline{\log_{k}(y)}$ and the eigenvectors can be chosen as the columns of the product of matrices

[TABLE]

Now we are ready to prove the following result.

Theorem 3.5**.**

Let $M$ be an invertible K80 Markov matrix with eigenvalues $1$, $x$, $y$, $y$ satisfying $x>0$ and $y\neq x$, and consider the family of real logarithms $Q(k,\alpha,\beta)$ for $M$ introduced in (10). If either $k\neq 0$ or $y<0$, an equality $Q(k,\alpha,\beta)=Q(k^{\prime},\alpha^{\prime},\beta^{\prime})$ implies $k=k^{\prime}$, $\alpha=\alpha^{\prime}$ and $\beta=\beta^{\prime}$. Otherwise, if $k=0$ and $y>0$ then $Q(0,\alpha,\beta)=\mathrm{Log}(M)$ for any $\alpha\in\mathbb{R},\beta\in\mathbb{R}_{>0}$.

Proof.

From Remark 3.4 we have that $(1,1,1,1)^{t}$, $(1,1,-1,-1)^{t}$ are eigenvectors with respective eigenvalues [math] and $\log(x)$. The other two eigenvalues may be complex:

[TABLE]

and the corresponding eigenvectors are

[TABLE]

Let us assume that $Q(k,\alpha,\beta)=Q(k^{\prime},\alpha^{\prime},\beta^{\prime})$ for some other choice of $k^{\prime},\alpha^{\prime},\beta^{\prime}$.

1st case: assume that the eigenvalues in (11) are complex. It is straightforward to check that this is the case if and only if $k\neq 0$ or if $y<0$. In any case, the four eigenvalues are different and so, the eigenspaces are 1-dimensional. From the equality of eigenvalues for $Q(k,\alpha,\beta)$ and $Q(k^{\prime},\alpha^{\prime},\beta^{\prime})$, we obtain that

[TABLE]

If $2\pi k+\theta_{y}=2\pi k^{\prime}+\theta_{y}$, we derive that $k=k^{\prime}$. Moreover, the corresponding eigenspaces must be equal, so $u_{+}(\alpha,\beta)=\lambda\;u_{+}(\alpha^{\prime},\beta^{\prime})$, for some $\lambda\in\mathbb{C}\setminus\{0\}$. From the first and third coordinates of these vectors we have $1+\alpha+\beta i=\lambda(1+\alpha^{\prime}-\beta^{\prime}i)$ and $1-\alpha-\beta i=\lambda(1-\alpha^{\prime}+\beta^{\prime}i)$ and by summing these two equations we derive that $\lambda=1$, and $\alpha=\alpha^{\prime}$ and $\beta=\beta^{\prime}$. If $2\pi k+\theta_{y}=-(2\pi k^{\prime}+\theta_{y})$, then $u_{+}(\alpha,\beta)=\lambda\;u_{-}(\alpha^{\prime},\beta^{\prime})$, for some $\lambda\in\mathbb{C}\setminus\{0\}$. Similarly as above, we deduce that $\lambda=1$. However, in this case, this implies that $\beta=-\beta^{\prime}$, which contradicts the assumption $\beta,\beta^{\prime}>0$.

2nd case: it remains to deal with the case where the eigenvalues in (11) are real, that is, $k=0$ and $y>0$. In this case, the eigenvalues of $Q(0,\alpha,\beta)$ are the principal logarithms of the eigenvalues of $M$. From the uniqueness of the principal logarithm, we infer that $Q(0,\alpha,\beta)=Log(M)$ for all $\alpha,\beta$.

∎

While the previous results list all the real logarithms with rows summing to zero of any K80 Markov matrix, we are mainly interested in those that are rate matrices (i.e. we need to restrict to non-negative values in the off-diagonal entries). The following results characterize the matrices $Q(k,\alpha,\beta)$ that are rate matrices.

Lemma 3.6**.**

With the notation of Theorem 3.5, we have that $Q(k,\alpha,\beta)$ is a rate matrix if and only if the following inequalities hold:

[TABLE]

Proof.

From Proposition 3.2 we have that $Q(k,\alpha,\beta)$ is a real matrix with rows summing to [math]. Thus we only need to characterize those $Q(k,\alpha,\beta)$ that have non-negative entries outside the diagonal. Let us define $\lambda=-\log(x)$ and $\mu=-\log|y|$. Note that $\lambda,\mu>0$ because $x,|y|\in(0,1]$ (see Remark 2.4). By computing $L_{0}$ as in Definition 3.1 we get that:

[TABLE]

On the other hand, given $A=\begin{pmatrix}1&0\\ \alpha&\beta\\ \end{pmatrix}\in GL_{2}(\mathbb{R})$ we obtain the following expression for $L_{1}(A)$:

[TABLE]

Recall that Q(k,\alpha,\beta)=L_{0}+\big{(}2\pi k+\theta_{y}\big{)}L_{1}(A). By looking at the off-diagonal entries, we get that $Q(k,\alpha,\beta)$ is a rate matrix if and only if:

[TABLE]

The first inequality above gives (12), while (13) follows by joining the second and third inequalities. ∎

In the following result, we prove that if $M$ is embeddable then $Q(0,0,1)$ is a Markov generator.

Theorem 3.7**.**

Let $M$ be an invertible K80 Markov matrix with eigenvalues $1$, $x$, $y$, $y$ satisfying $x>0$ and $y\neq x$. Then, if $Q(k,\alpha,\beta)$ is a rate matrix for some $k\in\mathbb{Z}$, $\alpha\in\mathbb{R}$ and $\beta\in\mathbb{R}_{>0}$ it holds that:

1. i)

$Q(l,\alpha,\beta)$* is a rate matrix for any integer $l\in I_{k}$ where*

[TABLE]

(we use the notation $\langle a,b\rangle$ to denote the closed interval delimited by $a$ and $b$, no matters if $a>b$ or $a<b$). 2. ii)

$Q(k,0,1)$* is a rate matrix.*

Proof.

We will prove that $Q(l,\alpha,\beta)$ and $Q(k,0,1)$ are Markov matrices by checking that they satisfy the inequalities (12) and (13) in Lemma 3.6.

1. i)

The proof is straightforward from Lemma 3.6 and the fact that $|2\pi l+\theta_{y}|\leq|2\pi k+\theta_{y}|$ for any $l\in I_{k}$ (note that $\theta_{y}$ is either 0 if $y>0$ or $\pi$ if $y<0$). 2. ii)

Since $Q(k,\alpha,\beta)$ is a rate matrix it follows from Lemma 3.6 that:

[TABLE]

Now, let us consider the real function $f(\beta)=\frac{1+\beta^{2}}{\beta}$ restricted to $\mathbb{R}_{>0}$. $f$ is continuous and a straightforward computation shows that $f$ has an absolute minimum at $\beta=1$. This concludes the proof.

∎

Now we are ready to prove the main result in this section, which characterizes embeddable K80 Markov matrices.

Corollary 3.8**.**

A K80 Markov matrix $M$ with eigenvalues $1$, $x$, $y$, $y$ is embeddable if and only if $Q(0,0,1)$ is a rate matrix, where

[TABLE]

In particular, if $y>0$ then $M$ is embeddable if and only if $\mathrm{Log}(M)$ is a rate matrix.

Proof.

If $x<0$ or $\det(M)=0$, then $M$ is not embeddable (see Remark 2.4). For $x=y>0$, $M$ is a JC69 matrix and it is known that a JC69 matrix is embeddable if and only if its eigenvalues are positive (Roca-Lacostena and Fernández-Sánchez,, 2018). In this case, $y>0$ and $Q(0,0,1)=\mathrm{Log}(M)$ is a rate matrix (see Theorem 2.5). Finally, if $x\neq y$ and $x>0$ the first claim follows from Proposition 3.2, Remark 3.3 and Theorem 3.7. To conclude the proof, note that if $y>0$ then $\theta_{y}=0$ thus $Q(0,0,1)=\mathrm{Log}(M)$ (see Definition 2.2). ∎

The following corollaries use Corollary 3.8 to characterize embeddable K80 Markov matrices in terms of its eigenvalues (Corollary 3.9) and in terms of its entries (Corollary 3.11).

Corollary 3.9**.**

Let $M$ be a K80 Markov matrix with eigenvalues $1$, $x$, $y$, $y$. Then:

1. (i)

If $y=0:$ $M$ is not embeddable. 2. (ii)

If $y>0:$ $M$ is embeddable if and only if $x\geq y^{2}$. In this case, $Log(M)$ is a rate matrix. 3. (iii)

If $y<0:$ $M$ is embeddable if and only if $e^{-2\pi}\geq x\geq y^{2}$.

Proof.

By Corollary 3.8 we know that $M$ is embeddable if and only if $Q(0,0,1)$ is a rate matrix. If $x\neq y$ it follows from Lemma 3.6 that $Q(0,0,1)$ is a rate matrix if and only if $\log(x)-2\log|y|\geq 0$ and $-\log(x)\geq 2\theta_{y}$. By computing the exponential of both sides of the expression $\log(x)\geq 2\log|y|$ we get $x\geq y^{2}$. Now, if $y<0$ we have $\theta_{y}=\pi$ and hence $-log(x)\geq 2\theta_{y}$ implies $e^{-2\pi}\geq x$. On the other hand, if $y>0$ we have that $\theta_{y}=0$ and hence $-\log(x)\geq 2\theta_{y}$ implies that $1\geq x$ . Note that this last constraint is redundant since the absolute value of the eigenvalues of any Markov matrix is bounded by 1 due to Perron-Frobenius theorem. ∎

Remark 3.10**.**

In the case of positive eigenvalues $x$ and $y$, we have that $M$ is embeddable if and only if its principal logarithm $\mathrm{Log}(M)$ is a rate matrix (Corollary 3.9). However, when the multiple eigenvalue $y$ is negative, $\mathrm{Log}(M)$ is never a rate matrix. In this case, the embeddability cannot be checked by looking at $\mathrm{Log}(M)$ and the previous corollary provides an embeddability criterion. In particular, the third item in the corollary above (or in Corollary 3.11 below) provides a non-zero measure subset formed by embeddable K80 matrices whose principal logarithm is not a rate matrix.

Corollary 3.11**.**

Let $M=K(1-b-2c,b,c,c)$. Then:

1. (i)

If ${2}c=1-2b$, $M$ is not embeddable. 2. (ii)

If ${2}c<1-2b$, $M$ is embeddable if and only if $c\leq\sqrt{b}-b\ (\leq\frac{1}{4})$. In this case, $Log(M)$ is a rate matrix. 3. (iii)

If ${2}c>1-2b$, $M$ is embeddable if and only if $\frac{1}{4}-\frac{e^{-2\pi}}{4}\leq c\leq\sqrt{b}-b\ (\leq\frac{1}{4})$.

Proof.

The claim follows directly from Corollary 3.9 by expressing the eigenvalues in terms of the entries using the bijection given in (3). ∎

Remark 3.12**.**

As claimed in the introduction, the results above together with Corollary 3.5 by Roca-Lacostena and Fernández-Sánchez, (2018) solve the embedding problem for any K81 matrix. Indeed, the embeddability of any K81 matrix with two repeated eigenvalues can be treated analogously to K80 matrices (one just needs to permute rows and columns accordingly). If the matrix has three repeated eigenvalues, then it is a JC69 matrix and the following argument solves the embeddability of JC69 matrices: there are no embeddable matrices with negative determinant (Culver,, 1966) and by Theorem 2.5, the principal logarithm of a JC69 matrix with positive determinant is always a rate matrix. Thus, a JC69 matrix is embeddable if and only if its determinant is positive.

4 Identifiability of rates for K80 Markov matrices

In this section we address the identifiability problem for K80 embeddable matrices. As a consequence of the results obtained in the previous section we provide a criterion to determine whether the rates of these matrices are identifiable or not. Furthermore, for those matrices with non-identifiable rates we determine how many Markov generators they admit.

Proposition 4.1**.**

Let $M$ be an embeddable K80 Markov matrix with eigenvalues $1$, $x$, $y$, $y$ and $x\neq y$. Then, the rates of $M$ are identifiable if and only if $Q(-1,0,1)$ is not a rate matrix.

Proof.

By Corollary 3.8 we know that $M$ is embeddable if and only if $Q(0,0,1)$ is a rate matrix.

Now assume that there are $\alpha\in\mathbb{R}$, $\beta\in\mathbb{R}_{>0}$ and $k\in\mathbb{Z}$ such that $Q(k,\alpha,\beta)$ is a Markov generator for $M$ different from $Q(0,0,1)$. In this case, Theorem 3.7 gives that $Q(-1,0,1)$ is also Markov generator because $-1$ belongs to the interval $I_{k}$, independently of the sign of $y$ (note that if $y>0$, the case $k=0$ is excluded by Theorem 3.5). Moreover, according to Theorem 3.5, $Q(0,0,1)$ and $Q(-1,0,1)$ are distinct Markov generators. ∎

Example 4.2**.**

In this example we show an embeddable K80 Markov matrix with positive eigenvalues and non-identifiable rates. Let us consider $M$ the K80 Markov matrix with eigenvalues $1$, $x=e^{-4\pi}$ and $y=e^{-2\pi}$ (with multiplicity 2). Rounding to the 10th decimal the entries of $M$ are:

[TABLE]

A straightforward computation shows that $\mathrm{Log}(M)$ is a Markov generator and hence $M$ is embeddable:

[TABLE]

Nonetheless, the rates of $M$ are not identifiable since there are other Markov generators for it:

[TABLE]

Note that this matrix has determinant equal to $e^{-8\pi}$ and is not close to the identity matrix $Id_{4}$. Indeed, in the Frobenius distance, $\|M-Id_{4}\|_{F}\approx 1.729893126$. At the other extreme, for any positive K80 matrix $K$, the limiting matrix $H:=\lim_{n\rightarrow\infty}K^{n}$ is the $4\times 4$ matrix whose entries are all $1/4$ and for the matrix above we have $\|M-H\|_{F}\approx 0.002640965$ (so $M$ is closer to $H$ than to $Id_{4}$).

Example 4.3**.**

Here we show an embeddable K80 Markov matrix with some negative eigenvalues and non-identifiable rates. Let us consider $M$ the K80 Markov matrix with eigenvalues $1$, $x=e^{-2\pi}$ and $y=-e^{-\pi}$ (with multiplicity 2). Rounding to the 10th decimal the entries of $M$ are:

[TABLE]

As we can see, $\mathrm{Log}(M)$ is not a real matrix:

[TABLE]

In spite of that, $Q(0,0,1)$ is a rate matrix, so $M$ is embeddable. Furthermore, $Q(-1,0,1)$ is also a Markov generator for $M$ and hence the rates are not identifiable:

$Q(0,0,1)=\footnotesize\begin{pmatrix}-\pi&0&0&\pi\\ 0&-\pi&\pi&0\\ \pi&0&-\pi&0\\ 0&\pi&0&-\pi\\ \end{pmatrix}$*   $Q(-1,0,1)=\footnotesize\begin{pmatrix}-\pi&0&\pi&0\\ 0&-\pi&0&\pi\\ {0}&{\pi}&-\pi&0\\ {\pi}&{0}&0&-\pi\\ \end{pmatrix}.$*

Unlike the previous example, in this case $Q(1,0,1)$ is not a rate matrix. Furthermore, $M$ is a K80 embeddable matrix with no K80 Markov generators (because $\mathrm{Log}(M)$ is not a rate matrix, see Theorem 2.5), which shows that rates do not necessarily satisfy the same symmetry constraints as probabilities do (Kimura,, 1980); see also (Roca-Lacostena and Fernández-Sánchez,, 2018).

With the notation of the previous example, note that this matrix has $\|M-Id_{4}\|_{F}\approx 1.781252133$, $\det(M)=e^{-4\pi}$, and $\|M-H\|_{F}\approx 0.06114223420$.

Remark 4.4**.**

In the examples 4.2 and 4.3, the Markov generators other than the principal logarithm are not K81 matrices, they belong to one of the Lie Markov models listed by Fernández-Sánchez et al., (2015), namely the model 3.3b. This is another 3-dimensional model, different from the K81 model, which contains the K80 model as well.

Theorem 4.5**.**

Let $M$ be an embeddable K80 Markov matrix with eigenvalues $1$, $x$, $y$, $y$ and $x\neq y$. Then the following holds:

1. a)

If $y>0$ and $x>e^{-4\pi}$ then $M$ has only one Markov generator, which is its principal logarithm. 2. b)

If $y>0$ and $x=e^{-4\pi}$ then $M$ has exactly $3$ generators: $Q(0,0,1)$ (which coincides with $\mathrm{Log}(M)$), $Q(1,0,1)$ and $Q(-1,0,1)$. 3. c)

If $y<0$ and $x=e^{-2\pi}$ then $M$ has exactly $2$ generators: $Q(0,0,1)$ and $Q(-1,0,1)$.

Otherwise, $M$ has infinitely many Markov generators.

Figure 1 illustrates Theorem 4.5 and Theorem 1.1.

Proof.

For ease of reading let us consider $\phi:\mathbb{R}\times\mathbb{R}_{>0}\rightarrow\mathbb{R}$, defined by \phi(\alpha,\beta)=\frac{\big{(}1+|\alpha|\big{)}^{2}+\beta^{2}}{\beta}. In the proof of Theorem 3.7 we already saw that $\phi$ has an absolute minimum at $\phi(0,1)=2$.

Note that, since $M$ is embeddable, we have $x\geq y^{2}$ (Corollary 3.9) and hence inequality (12) is satisfied for any $\alpha,\beta$ such that $\alpha^{2}+\beta^{2}=1$ or for $k=0$ when $y>0$.

Let $Q(k,\alpha,\beta)$ be a Markov generator for $M$.

1. a)

We know that $Q(-1,0,1)$ is not a rate matrix by Lemma 3.6 (it does not satisfy (13)). Hence, by Proposition 4.1 the rates of $M$ are identifiable. Furthermore, due to Corollary 3.8 we have that the only Markov generator of $M$ must be its principal logarithm. 2. b)

Since $y>0$ we have that $\theta_{y}=0$. Hence, it follows from inequality (13) in Lemma 3.6 that $4\pi\geq 2|k|\pi\;\phi(\alpha,\beta)$. Using that $\phi(\alpha,\beta)>2$ for $(\alpha,\beta)\neq(0,1)$ we get that the inequality holds if and only if $k=0$ (and hence $Q(0,\alpha,\beta)=\mathrm{Log}(M)$) or $\alpha=0$, $\beta=1$ and $|k|=1$. Note that all these solutions do also satisfy inequality (12). 3. c)

Since $y<0$ we have that $\theta_{y}=\pi$. Hence, it follows from inequality (13) in Lemma 3.6 that 2\pi\geq\big{|}(2k+1)\pi\;\phi(\alpha,\beta)\big{|}. Using that $\phi(\alpha,\beta)>2$ for $(\alpha,\beta)\neq(0,1)$ we get that the inequality holds if and only if $k=0$ (and hence $Q(0,\alpha,\beta)=\mathrm{Log}(M)$) or $\alpha=0$, $\beta=1$ and $|2k+1|\leq 1$. Note that all these solutions do also satisfy inequality (12).

If $M$ does not lie in the already covered cases, note that if $y>0$ we have $x<e^{-4\pi}$. In particular, $-\log(x)>4\pi=2\pi|-1|\;\phi(0,1)$ and inequality (13) is satisfied for $Q(-1,0,1)$. Furthermore, inequality (13) is satisfied for any $(\alpha,\beta)$ close enough to $(0,1)$. Hence, it follows from Lemma 3.6 that $Q(-1,\alpha,\beta)$ is a Markov generator of $M$ for any $(\alpha,\beta)$ close enough to $(0,1)$ such that $\alpha^{2}+\beta^{2}=1$. The same argument does also work for $y<0$. ∎

Remark 4.6**.**

The determinant of $M$ is related to the expected number of nucleotide substitutions in the Markov process ruled by $M$. In a phylogenetic tree, the length of a branch representing an evolutionary process between an ancestral species and a descendant species is usually measured as the expected number of nucleotide substitutions per site. If this process is ruled by a substitution matrix with uniform stationary distribution then this expected number of substitutions can be approximated by $l(M):=-\frac{1}{4}log(\det(M))$ (for a precise formulation see Barry and Hartigan,, 1987). For example, for the identity matrix one has $l(Id_{4})=0$ and for the limiting matrix $H$ introduced in Example 4.2, $l(H)=\infty$. We come back to this biological concept in the Discussion but for the moment we give conditions on the identifiability in terms of $l(M)$.

Corollary 4.7**.**

Let $M$ be an embeddable K80 Markov matrix with eigenvalues $1$, $x$, $y$, $y$ with $y\neq x>0$. Then,

1. (a)

If $y>0$ and $l(M)<{2\pi}$, the rates of $M$ are identifiable. 2. (b)

If $y<0$ and $l(M)<{\pi}$, the rates of $M$ are identifiable.

Moreover, these bounds are sharp (tight examples are given in 4.2 with $l(M)=2\pi$ and in 4.3 with $l(M)=\pi$).

Proof.

Theorem 4.5 shows that the eigenvalues $x,y$ of embeddable K80 matrices must satisfy $x\geq y^{2}$ and hence their determinant $xy^{2}$ must be smaller than or equal to $x^{2}$. On the other hand, Theorem 4.5 provides a bound on $x$ (depending on the sign of $y$) to determine whether an embeddable matrix has identifiable rates or not. More precisely, K80 embeddable matrices with non-identifiable rates have determinant at most $e^{-8\pi}$ if all the eigenvalues are positive or determinant at most $e^{-4\pi}$ if they have a repeated negative eigenvalue. Equivalently, their branch length is greater than or equal to $2\pi$ or $\pi$, respectively.

As we noted in the Examples 4.2 and 4.3 above, the K80 embeddable matrices provided there have non-identifiable rates and determinant $e^{-8\pi}$ and $e^{-4\pi}$, respectively. ∎

In the following corollary we restate Theorem 4.5 in terms of the entries of the Markov matrix.

Corollary 4.8**.**

Let $M=K(1-b-2c,b,c,c)$ with $b\neq c$ be an embeddable K80 Markov matrix. Then the following holds:

1. i)

If $2c<1-2b$ and $c<\frac{1}{4}-\frac{e^{-4\pi}}{4}$ then $M$ has only one Markov generator, which is its principal logarithm. 2. ii)

If $2c<1-2b$ and $c=\frac{1}{4}-\frac{e^{-4\pi}}{4}$ then $M$ has exactly $3$ generators: $Q(0,0,1)$ (which coincides with $\mathrm{Log}(M)$), $Q(1,0,1)$ and $Q(-1,0,1)$. 3. iii)

If $2c>1-2b$ and $c=\frac{1}{4}-\frac{e^{-2\pi}}{4}$ then $M$ has exactly $2$ generators: $Q(0,0,1)$ and $Q(-1,0,1)$.

Otherwise, $M$ has infinitely many Markov generators.

Proof.

The claim follows directly from Theorem 4.5 by expressing the eigenvalues in terms of the entries using the bijection given in (3). ∎

We can now provide a proof for Theorem 1.1:

Proof of Theorem 1.1. Theorem 1.1 is the summary of the results in Corollary 3.11 and Corollary 4.8. $\square$

We turn our attention to matrices whose diagonal entries are the largest entries in each column. These are called diagonal largest in column matrices, briefly DLC, and are related to matrix parameter identifiability in phylogenetics (see Chang,, 1996). A straightforward computation shows that

[TABLE]

Note that the Markov matrix of Example 4.2 is DLC, thus there are DLC matrices whose rates are not identifiable.

Corollary 4.9**.**

The following holds:

1. (a)

A K80 DLC matrix is embeddable if and only if its principal logarithm is a rate matrix. 2. (b)

There is an open set of embeddable DLC matrices whose rates are not identifiable.

Proof.

It follows from the first inequality in (14) that DLC K80 Markov matrices lie in case (b) of Theorem 1.1. The second statement follows from case (b) iii) in the same Theorem.

5 Measuring the set of embeddable and rate identifiable matrices

In this section we will study how many K80 Markov matrices are embeddable and also how many of those matrices have identifiable rates. We will solve these questions when restricted to some well defined subsets of K80 matrices. We will proceed by using a parametrization of all K80 matrices and computing the volumes (actually areas) needed in that set. We will consider the following regions:

1. i)

$\Delta$: The set of all $K80$ Markov matrices. 2. ii)

$\Delta_{\rm{DLC}}$: The set of all DLC matrices in $\Delta$ (see Remark 4.9). 3. iii)

$\Delta_{+}$: The set of all $M\in\Delta$ in the connected component of the identity with respect to vanishing determinant (that is, all K80 Markov matrices with only positive eigenvalues). 4. iv)

$\Delta_{\rm{dd}}$: The set of all $M\in\Delta$ such that the probability of not mutating is higher than the probability of mutating, i.e. those matrices such that $a>b+2c$. These matrices are said to be diagonally-dominant and it is well known that they have identifiable rates (if embeddable) (see Cuthbert,, 1972).

In the previous sections we have been using two different parameterizations of the model (see Equation (3) and Figure 1), one in terms of the eigenvalues of the Markov matrix ($x$ and $y$) and the other in terms of its entries ($b$ and $c$) . The first parametrization can be used to easily describe $\Delta_{+}$ whereas the second provides an easier description of $\Delta_{\rm{DLC}}$ and $\Delta_{\rm{dd}}$. Using these two parametrizations, we obtain that $\Delta\supset\Delta_{\rm{DLC}}\supset\Delta_{+}\supset\Delta_{\rm{dd}}$. Figures 2 and 3 illustrate these inclusions in terms of the eigenvalues and the entries of the Markov matrix, respectively.

Note that, since the bijection $\varphi$ in (3) is a linear map, the relative volumes of embeddable matrices will not depend on the chosen parameters (entries or eigenvalues). Actually, $\det(D(\varphi))=8$, so the volume of any subset will be eight times greater in the eigenvalues’ parametrization. We decided to use the parametrization in terms of the eigenvalues to compute the volumes we want to know because the expressions appearing in the computations are simpler.

For a clearer picture of embeddability and rate identifiability of those subsets of K80 matrices, intersect Figure 1 and Figures 2 and 3.

Proposition 5.1**.**

Let $\Delta_{\rm{emb}}$ denote the set of embeddable K80 Markov matrices, and let $\Delta_{\rm{idf}}\subset\Delta_{\rm{emb}}$ be the subset of matrices with identifiable rates. Using the parametrization of K80 Markov matrices in terms of its eigenvalues, the following holds: a) $V(\Delta)=2$; b) $V(\Delta_{\rm{DLC}})=10/12$; c) $V(\Delta_{+})=3/4$; d) $V(\Delta_{\rm{dd}})=1/2$; e) $V(\Delta_{\rm{emb}})=\frac{2(1+e^{-3\pi})}{3}$; f) $V(\Delta_{\rm{idf}})=\frac{2(1-e^{-6\pi})}{3}$; g) $V(\Delta_{\rm{emb}}\cap\Delta_{\rm{DLC}})=V(\Delta_{\rm{emb}}\cap\Delta_{+})=2/3$; h) $V(\Delta_{\rm{emb}}\cap\Delta_{\rm{dd}})=V(\Delta_{\rm{idf}}\cap\Delta_{\rm{dd}})=\frac{7-4\sqrt{2}}{3}$.

Proof.

1. a)

$\Delta$ is the triangle with vertices $(-1,0),\ (1,1)\text{ and }(1,-1)$ which has area 2. 2. b)

$\Delta_{\rm{DLC}}$ is the polygon with vertices $(-1/3,1/3),\ (1,1),\ (1,0)\text{ and }(0,0)$ which has area 10/12. 3. c)

$\Delta_{+}$ is the trapezoid with vertices $(0,1/2),\ (1,1),\ (1,0)\text{ and }(0,0)$ which has area 3/4. 4. d)

$\Delta_{\rm{dd}}$ is the triangle with vertices $(0,1/2),\ (1,1)\text{ and }(1,0)$ which has area 1/2. 5. e)

It follows from Corollary 3.9 that $V(\Delta_{\rm{emb}})=\displaystyle\int_{0}^{1}\int^{1}_{y^{2}}1\ dxdy\ +\ \ \int_{-e^{-\pi}}^{0}\int_{y^{2}}^{e^{-2\pi}}1\ dxdy=\frac{2(1+e^{-3\pi})}{3}$. 6. f)

It follows from Theorem 4.5 that $V(\Delta_{\rm{idf}})=\displaystyle\int_{0}^{1}\int^{1}_{y^{2}}1\ dxdy\ -\ \ {\int_{0}^{e^{-2\pi}}}\int_{y^{2}}^{e^{-4\pi}}1\ dxdy=\frac{2(1-e^{-6\pi})}{3}$. 7. g)

As shown in the proof of Corollary 4.9, the entries of a DLC K80 Markov matrix must satisfy $1-b-2c>b$ and $1-b-2c>c$. Hence, the repeated eigenvalue of a K80 DLC Markov matrix $y=1-2b-2c$ (see equation (2)), must be positive. Since embeddable matrices have positive determinant, we have that $\Delta_{\rm{emb}}\cap\Delta_{\rm{DLC}}=\Delta_{\rm{emb}}\cap\Delta_{+}$ (see Figure 2). It follows from Corollary 3.9 that the volume of this set can be computed as $\displaystyle\int_{0}^{1}\sqrt{x}\ dx=2/3$. 8. h)

The set of diagonally-dominant matrices is the triangle with vertices $(0,0.5),\ (1,1)$ and $(1,0)$. It is known that diagonally-dominant Markov matrices have only one real logarithm (Cuthbert,, 1972) and hence the first equality follows. Furthermore the points where the curve $y^{2}=x$ intersects with the boundary of $\Delta_{\rm{dd}}$ are $(3-2\sqrt{2},-1+\sqrt{2})$ and $(1,1)$, thus it follows from Corollary 3.9 that the volume can be computed as$\displaystyle\int^{1}_{\sqrt{2}-1}\int^{1}_{y^{2}}1\ dxdy+2(\sqrt{2}-1)^{2}=\frac{7-4\sqrt{2}}{3}$

∎

The table 1 shows some relative volumes of the regions defined at the beginning of this section, and is included here for quick reference. The computations involved are straightforward from the values of the volumes of the preceding proposition.

6 Discussion

We have studied the embeddability and the identifiability of mutation rates for the Markov matrices in the K80 model of nucleotide substitution. With the results of the present paper and (Roca-Lacostena and Fernández-Sánchez,, 2018), the problem of embeddability is completely solved for the K81 model and its submodels K80 and JC69. One of the relevant results of the present paper is proving the existence of an open subset of embeddable K80 matrices whose principal logarithm is not a rate matrix (see Remark 3.10). Moreover, we have also provided an open set of matrices in the model that have infinitely many Markov generators (Theorem 1.1 $b.iii)$ and $c.ii)$). Note that this might lead to confusing results when trying to infer parameters as it is usually done in phylogenetics via a maximum likelihood approach. We further develop this issue in what follows.

In phylogenetic trees that evolve under a Markov evolutionary model one usually assigns a branch length to the edges of the tree (see Remark 4.6). This length accounts for the expected number of substitutions per site that have occurred along that edge. It is well known that, if a K80 Markov matrix $M$ has governed evolution on an edge, then the expected number of elapsed substitutions per site can be approximated by $l(M)=\frac{-1}{4}\log(\det(M))$ as long as $M$ is a product of Markov matrices close to the identity (Barry and Hartigan,, 1987). According to Corollary 4.7, there are embeddable K80 matrices $M$ with $l(M)=\pi$ whose rates are not identifiable and whose principal logarithm is not a rate matrix (Example 4.3). Therefore, there are K80 matrices with branch length $\pi$ whose embeddability property cannot be decided by looking at the principal logarithm (this is a usual practice when deciding embeddability, see Verbyla et al., (2013), and might lead to erroneous results as shown here). The expected number of substitutions per site is related to molecular rate: different species evolve at different molecular rates and therefore, whether an expected number of substitutions equal to $\pi$ is large or not, depends on the species under consideration. This notion can be also related to (astronomical) time if a molecular clock can be assumed. As an example, the molecular rate is estimated to be in the range of 0.1-10 substitution per site per million year for bacteria (Duchene et al.,, 2016), 0.03-0.1 for birds or 0.2-0.8 for humans (Ho et al.,, 2007). Therefore, a branch of length $\pi$ could represent from $\pi/10$ to $\pi/0.03\approx 104.72$ million years, depending on the species considered (at least for the examples provided above).

On the other hand, as noted in Corollary 3.8, these cases where embeddability cannot be decided by looking at the principal logarithm correspond to negative multiple eigenvalue $y$ (Remark 3.10). Although these cases might not seem biologically realistic, there are embeddable K80 matrices with negative eigenvalues that are closer to the identity matrix than some other embeddable K80 matrices with positive eigenvalues. For more general models, one might be able to find examples of embeddable matrices even closer to the identity matrix with such a behaviour (work in progress).

We have obtained that K80 DLC matrices are embeddable if and only if its principal logarithm is a rate matrix (Remark 4.9). In parameter estimation in phylogenetics, it might be relevant to restrict to a subset of matrices where identifiability of the substitution parameters is guaranteed (for a discussion see Zou et al.,, 2011 and Kaehler et al.,, 2015). The set of all DLC matrices is one of these subsets (see Chang,, 1996). We have found an open subset of DLC matrices whose rates are not identifiable (see Corollary 4.9), which actually have infinitely many Markov generators. This implies that different continuous time processes could lead to the same final observations, so one has to be careful when doing maximum likelihood estimation. The DLC matrix that is closest to the identity and has non-identifiable rates corresponds to branch length $2\pi$, see Example 4.2 and Theorem 4.5.

In a work in progress, we are studying a more general Markov model which allows complex eigenvalues and therefore opens the door to exploring ill-behavior embeddable matrices (in terms of principal log not being real) that are closer to the identity matrix. Moreover, initial studies for the general Markov models suggest the existence of an open subset of embeddable Markov matrices with non-identifiable rates. We will address this case in a forthcoming paper. The identifiability of rates for JC69 matrices is a problem not solved in the present paper (the case of an eigenvalue with multiplicity 3 adds some technical difficulties), but we expect to solve it in future work.

The computation of relative volumes carried out in Section 5 shows that less than 35% of the Markov matrices within the K80 model are embeddable. This suggests that restricting to the continuous-time approach within this model supposes a strong restriction for the inference of the Markov matrix that rules a given evolutionary process. At the same time, this restriction guarantees a biological realism in the modeling process (see Theorem 2.2 by Roca-Lacostena and Fernández-Sánchez, (2018)) and might be an advantage for reconstruction methods that estimate the parameters of the model. On the other hand, restricting to continuous-time models might lead to overestimation of genetic distances for nonstationary data (Kaehler et al.,, 2015) and, in general, the choice of an accurate model for the data is indispensable. The computations in Section 5 also exhibit that even if K80 embeddable matrices with non-identifiable rates describe a subset of positive measure, this case is quite marginal and represents less than 1 over 10.000 embeddable matrices.

Author contribution

MC and JFS conceived the project, revised the proofs and computations and drafted part of the manuscript. JRL wrote the core of the manuscript and worked out the proofs and computations. All authors read, revised and approved the final manuscript.

Acknowledgements

All authors are partially funded by AGAUR Project 2017 SGR-932 and MINECO/FEDER Projects MTM2015-69135 and MDM-2014-0445. J Roca-Lacostena has received also funding from Secretaria d’Universitats i Recerca de la Generalitat de Catalunya (AGAUR 2018FI_B_00947) and European Social Funds. The authors would like to express their gratitude to Jeremy Sumner for his remarks and interesting conversations on the topic.