Models of sequence evolution 3

I wrote a couple of long posts about models of sequence evolution (here and here), in which (and it may not be obvious) we solved the problem of how to obtain:

P = eQt

for any t, by getting eigenvalues and eigenectors, using both R and Python. That's outstanding.

If you look carefully, you might notice something. There are just two terms in Q for the JC69 model:

-0.01 (on diagonal) and 0.003333 (off diagonal)

P has only two terms as well and they turn out to be (for t = 0.1):

0.999 and 0.000333

That is, the off-diagonal terms are apparently just multiplied by t, and naturally the diagonal is 1 - 3 times the other term. There was really no need for all the exponentiation of matrices, etc.

How does this happen? The formulas for the two terms in JC69 are:

PNN(t) = 1/4 + 3/4 exp(-4λt) PNM(t) = 1/4 - 1/4 exp(-4λt)

But going back to the infinite series expansion for e (or just remembering the approximation for small x):

ex = 1 + x1/1! + x2/2! + .. ex = 1 + x

So the exponential is just 1 - 4λt and we obtain:

PNN(t) = 1/4 + 3/4 (1 - 4λt) = 1 - 3λt PNM(t) = 1/4 - 1/4 (1 - 4λt) = λt

Let's take another look at the formulas for these terms in the TN93 model. Look at only the first row of P. We first consider the two exponential terms:

e2 = exp(-βt) e4 = exp[-(πYα1 + πRβ)t]

Approximating for small t:

e2 = 1 - βt e4 = 1 - (πYα1 + πRβ)t

Look what happens to P13 and P14:

P13 = &piA(1 - e2) = &piAβt

This is t times Q13. Exactly the same thing happens for P14. P12 is a bit more complicated, but it's just algebra:

P12 = πC + (πCπR/πY)e2 - (πC/πY)e4

The second term is:

πCπR/πY (1 - βt) πCπR/πY - (πCπR/πY)βt

The third term is:

πC/πY (-1 + πYα1t + πRβt)

We can see that (considering both the second and the third term) there are two copies of: (πCπR/πY)βt with opposite sign, leaving:

P12 = πC + πCπR/πY - πC/πY + πCα1t

The last term is t times Q12. So I'm betting that the rest will go away. We have:

πC + πCπR/πY - πC/πY

but

πY + πR = 1 πR - 1 = -πY

so

πC + πCπR/πY - πC/πY = πC + πC/πY(πR - 1) = πC - πC = 0

P12 = πCα1t = t Q12

How about that!