Jukes-Cantor (7)

I wrote a little script that derives the eigenvalues and eigenvectors for a simple Q matrix, and then tests whether U Λ U-1 = Q, as we said last time.

It looks good, except that I haven't been able to actually do the exponentiation part yet.

import numpy as np import math a = 0.00333333 Q = np.array( [[ -3*a, a, a, a], [ a, -3*a, a, a], [ a, a, -3*a, a], [ a, a, a, -3*a]]) evals, evecs = np.linalg.eig(Q) U = evecs L = np.diag(evals) Ui = np.linalg.inv(U) temp = np.dot(U,L) P = np.dot(temp,Ui) names = ['Q','evals','evecs','U','L','Ui','P'] for i,v in enumerate([Q,evals,evecs,U,L,Ui,P]): print names[i] for row in v: try: temp = [str(round(e,3)).rjust(5) for e in row] print ' '.join(temp) except: print row print

Q -0.01 0.003 0.003 0.003 0.003 -0.01 0.003 0.003 0.003 0.003 -0.01 0.003 0.003 0.003 0.003 -0.01 evals -0.01333332 -1.30104260698e-18 -0.01333332 -0.01333332 evecs -0.866 0.5 0.269 0.024 0.289 0.5 0.624 0.314 0.289 0.5 -0.182 -0.818 0.289 0.5 -0.711 0.481 U -0.866 0.5 0.269 0.024 0.289 0.5 0.624 0.314 0.289 0.5 -0.182 -0.818 0.289 0.5 -0.711 0.481 L -0.013 0.0 0.0 0.0 0.0 -0.0 0.0 0.0 0.0 0.0 -0.013 0.0 0.0 0.0 0.0 -0.013 Ui -0.866 0.542 0.235 0.089 0.5 0.5 0.5 0.5 0.0 0.789 -0.102 -0.688 0.0 0.322 -0.811 0.489 P -0.01 0.003 0.003 0.003 0.003 -0.01 0.003 0.003 0.003 0.003 -0.01 0.003 0.003 0.003 0.003 -0.01