Building a Phylogenetic Tree From Sequence Data

Distance-based tree building is three steps: compute pairwise distances, correct them, and join neighbours greedily. The second step is the one that is usually skipped and the one that decides whether the tree means anything.

From sequences to distances

Everything starts from a multiple sequence alignment, and the alignment is a modelling step, not a preprocessing step. Which columns are homologous is an inference, and a tree built on a bad alignment is a tree of the alignment’s errors. Columns containing gaps are typically either removed entirely or handled explicitly, and that choice alone can move branches.

Given the alignment, the simplest distance between two sequences is the p-distance: the fraction of aligned positions at which they differ. Count 12 differences over 200 comparable positions and p is 0.06. It is easy to compute and it is systematically wrong, in a direction that gets worse the more divergent the sequences are.

Alignment-free distances exist and avoid the alignment step entirely by comparing k-mer profiles between sequences — the same representation used for k-mer based sequence classification. They scale to whole genomes, where a multiple alignment is impractical, and they are the standard route for genome-scale phylogenies. What they give up is the site-level correspondence that makes the correction below possible: a k-mer distance has no notion of substitutions per site, so it cannot be corrected for multiple hits in the way described next, and its relationship to evolutionary time is empirical rather than modelled.

Why raw differences are not distances

The problem is multiple hits. A site that changed from A to G and later from G to A shows no difference today, and a site that changed twice to different bases shows one difference where two substitutions occurred. Observed differences therefore undercount actual substitutions, and the undercount compounds: as sequences diverge, p saturates toward 0.75 for DNA — the value expected for two random sequences over four equal-frequency bases — while the true number of substitutions keeps rising without bound.

A substitution model corrects for this. The simplest, from Jukes and Cantor in 1969, assumes all substitutions are equally likely and all bases equally frequent, and inverts the saturation analytically:

d = -(3/4) * ln(1 - (4/3) * p)

  p = 0.05  ->  d = 0.0518   ( +4% )
  p = 0.20  ->  d = 0.2326   ( +16% )
  p = 0.50  ->  d = 0.8240   ( +65% )
  p = 0.70  ->  d = 2.0313   ( +190% )
  p -> 0.75 ->  d -> infinity

The correction is negligible for close sequences and enormous for distant ones, which is exactly why an uncorrected tree looks fine on a cluster of related strains and collapses when you add an outgroup. Richer models — allowing transitions and transversions to differ, allowing unequal base frequencies, allowing rate variation across sites with a gamma distribution — correct further, and the divergence at which they matter is lower than most people expect.

The property all this is chasing is additivity: distances are usable by a tree-building algorithm when the distance between two taxa equals the sum of the branch lengths on the path between them. Corrected distances are approximately additive; raw p-distances are not, and neighbour joining assumes they are.

The neighbour-joining criterion

Neighbour joining, published by Naruya Saitou and Masatoshi Nei in Molecular Biology and Evolution in 1987, repeatedly joins the pair of taxa that most reduces the total length of the tree. The subtlety, and the reason it beats naive clustering, is that it does not join the closest pair.

Joining the closest pair — which is what a simple agglomerative clustering does — assumes a molecular clock, that all lineages evolve at the same rate. When they do not, two slowly evolving taxa can be closer to each other than either is to its true sibling, and clustering groups them wrongly. Neighbour joining corrects each pairwise distance by how far both members are from everything else:

S_i = sum over all k of d(i,k)          // divergence of i

Q(i,j) = (n - 2) * d(i,j) - S_i - S_j

join the pair minimising Q, then:

  delta(i,u) = 0.5 * d(i,j) + (S_i - S_j) / (2*(n-2))
  delta(j,u) = d(i,j) - delta(i,u)
  d(u,k)     = 0.5 * (d(i,k) + d(j,k) - d(i,j))

Subtracting the two divergence sums is the whole trick. A taxon on a long branch is far from everything, so its large S penalises every pair it appears in, and it stops being joined to whatever it happens to be least far from.

Four taxa, worked

Take four taxa with this corrected distance matrix, n equal to 4:

        A      B      C      D          S_i (row sums)
  A     -    0.10   0.40   0.45        0.95
  B   0.10     -    0.40   0.45        0.95
  C   0.40   0.40     -    0.15        0.95
  D   0.45   0.45   0.15     -         1.05

Q(A,B) = 2(0.10) - 0.95 - 0.95 = 0.20 - 1.90 = -1.70
Q(C,D) = 2(0.15) - 0.95 - 1.05 = 0.30 - 2.00 = -1.70
Q(A,C) = 2(0.40) - 0.95 - 0.95 = 0.80 - 1.90 = -1.10
Q(A,D) = 2(0.45) - 0.95 - 1.05 = 0.90 - 2.00 = -1.10
Q(B,C) = -1.10        Q(B,D) = -1.10

A-B and C-D tie at minus 1.70, which is not a defect: the matrix is symmetric between the two pairs, so the tie is the algorithm correctly reporting that both joins are equally supported. Break it arbitrarily and join A and B at a new node u:

delta(A,u) = 0.5(0.10) + (0.95 - 0.95)/(2*2) = 0.05
delta(B,u) = 0.10 - 0.05                       = 0.05

d(u,C) = 0.5(0.40 + 0.40 - 0.10) = 0.35
d(u,D) = 0.5(0.45 + 0.45 - 0.10) = 0.40

remaining matrix on {u, C, D}:
        u      C      D
  u     -    0.35   0.40
  C   0.35     -    0.15
  D   0.40   0.15     -

three taxa resolve directly to a star with centre v:
  delta(C,v) = 0.5(0.15 + 0.35 - 0.40) = 0.05
  delta(D,v) = 0.5(0.15 + 0.40 - 0.35) = 0.10
  delta(u,v) = 0.5(0.35 + 0.40 - 0.15) = 0.30

tree:  ((A:0.05, B:0.05):0.30, C:0.05, D:0.10)

Read the result. A and B are a tight pair on short branches; C and D are a pair separated by a long internal branch of 0.30; and D sits on a branch twice as long as C, which is the algorithm recovering that D was slightly further from everything in the original matrix. Notice that the output is unrooted — the tree describes relationships, not direction — and rooting requires either an outgroup known to be outside the group of interest, or a molecular-clock assumption you have to justify.

Where distance methods stop

The step from alignment to distance matrix discards almost all the information. A single number per pair cannot express that two taxa share a rare substitution at one particular site, which is exactly the evidence that resolves difficult relationships. That is why distance-based trees are fast, deterministic and used as starting points rather than as final answers.

• Maximum likelihood and Bayesian methods score whole trees. They compute, under an explicit substitution model, the likelihood of the alignment given a tree with branch lengths, and search tree space for the best. Every column contributes. The cost is that tree space grows super-exponentially with taxon count, so the search is heuristic and can get stuck.
• A tree without support values is not a result. Bootstrap resampling rebuilds the tree from alignment columns sampled with replacement, many times, and reports how often each branch recurs. Branches appearing in a small minority of replicates should be read as unresolved, and neighbour joining will always give you a fully bifurcating tree regardless of whether the data supports one.
• Long-branch attraction. Two rapidly evolving lineages accumulate so many substitutions that they share changes by chance and group together in the tree despite being unrelated. Distance correction reduces it and does not eliminate it; the diagnostics are adding taxa that break up the long branches and checking whether the grouping survives a model that allows rate variation across sites.
• A gene tree is not a species tree. Individual genes have their own histories through recombination, incomplete lineage sorting and horizontal transfer. Trees from different genes in the same organisms genuinely disagree, and reconciling them is a separate problem rather than an error in any one tree.