Solution (source code)

= Solution

Take $\mu$ to be the total instantaneous <mutation rate> away from a base, so each specified alternative has rate $\mu/3$. Under the <Jukes-Cantor model>, the uniform component is stationary and the three contrast components decay as $e^{-4\mu t/3}$, giving
$$
P_{ab}(t)=\begin{cases}
\frac14+\frac34e^{-4\mu t/3},&a=b,\\
\frac14(1-e^{-4\mu t/3}),&a\ne b.
\end{cases}
$$
An implementation is:

``
def jukes(a, b, t, mu):
    decay = np.exp(-4*mu*t/3)
    return .25 + .75*decay if a == b else .25*(1-decay)

print(f"{jukes(2,2,1,.5):.3f} {jukes(2,3,1,.5):.3f}")
``

\b[The two probabilities are $0.635$ and $0.122$.] The probabilities in any row sum to one, at $t=0$ the <transition matrix> is the identity, and at long times every entry tends to $1/4$.

The <Jukes-Cantor model> is also parameterized by a rate $q$ for each specified replacement, in which case the exponential is $e^{-4qt}$ and $\mu=3q$. The total-rate convention is used consistently below; it yields an optimum within the stated search interval, whereas the corresponding per-replacement rate is about $0.0522$.