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:
But going back to the infinite series expansion for e (or just remembering the approximation for small x):
So the exponential is just
1 - 4λtand we obtain:
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:
Approximating for small t:
Look what happens to
This is t times
Q13. Exactly the same thing happens for
P12is a bit more complicated, but it's just algebra:
The second term is:
The third term is:
We can see that (considering both the second and the third term) there are two copies of:
(πCπR/πY)βtwith opposite sign, leaving:
The last term is t times
Q12. So I'm betting that the rest will go away. We have:
How about that!