DistanceCalculator identity model inflates the distance when skip_letters is used #5332

Setup

I am reporting a problem with Biopython version, Python version, and operating
system as follows:

import sys; print(sys.version)
import platform; print(platform.python_implementation()); print(platform.platform())
import Bio; print(Bio.__version__)

3.12.10 (tags/v3.12.10:0cc8128, Apr  8 2025, 12:21:36) [MSC v.1943 64 bit (AMD64)]
CPython
Windows-11-10.0.26200-SP0
1.89.dev0

This also affects the released Biopython 1.88 (the same code is present in the
biopython-188 tag). numpy 2.x is required for DistanceCalculator to load a
substitution matrix, but the bug below reproduces with model="identity",
which needs no matrix.

Expected behaviour

DistanceCalculator._pairwise() documents that it returns "a value between 0
(identical sequences) and 1". Two identical sequences that contain gap columns
should therefore give a distance of 0.0, and the number of gap columns should
not affect the result at all, because gaps are explicitly excluded via
skip_letters.

Actual behaviour

With model="identity" the numerator and the denominator are built from
different sets of columns.

Bio/Phylo/TreeConstruction.py, DistanceCalculator._pairwise():

if self.scoring_matrix is None:
    # Score by character identity, not skipping any special letters
    score = sum(
        l1 == l2
        for l1, l2 in zip(seq1, seq2)
        if l1 not in self.skip_letters and l2 not in self.skip_letters
    )
    max_score = len(seq1)          # <-- full length, gaps included

score skips the positions named in skip_letters, but max_score is the full
sequence length, so every skipped column still inflates the denominator and the
returned distance is scaled by len(seq1) rather than by the number of positions
that were actually scored.

The scoring-matrix branch of the very same function does not have this problem,
because it accumulates max_score only over non-skipped positions. The two
branches therefore disagree on identical input.

This affects model="identity" combined with an explicit skip_letters
argument, which is the documented backward-compatibility path (see the
# Shim for backward compatibility (#491) in DistanceCalculator.__init__).
The default DistanceCalculator() / model="identity" sets skip_letters = (),
so the default path is unaffected and existing tests do not cover this.

Steps to reproduce

from io import StringIO

from Bio import Align, AlignIO
from Bio.Phylo.TreeConstruction import DistanceCalculator

# 1. Two IDENTICAL sequences, two of the six columns skipped.
aln = AlignIO.read(StringIO(">Alpha\nAC--GT\n>Gamma\nAC--GT"), "fasta")
dmat = DistanceCalculator("identity", skip_letters=("-", "*")).get_distance(aln)
print(repr(dmat["Alpha", "Gamma"]))   # actual:   0.3333333333333333
                                      # expected: 0.0

The four non-skipped columns all match, so the score is 4 and the maximum score
should be 4, giving 1 - 4/4 == 0.0. Instead max_score is len(seq1) == 6,
giving 1 - 4/6 == 0.3333....

The two branches of the same function disagree on the same data:

aln = AlignIO.read(StringIO(">Alpha\nAC--GTACGT\n>Gamma\nAC--GTACGT"), "fasta")
ident = DistanceCalculator("identity", skip_letters=("-", "*")).get_distance(aln)
# 'trans' takes the scoring-matrix branch, which excludes gaps from max_score
matrix = DistanceCalculator("trans", skip_letters=("-", "*")).get_distance(aln)
print(repr(ident["Alpha", "Gamma"]), repr(matrix["Alpha", "Gamma"]))
# actual:   0.2 0.0        <- the two branches of _pairwise() disagree
# expected: 0.0 0.0

And the distance ends up measuring gap columns rather than mismatches:

one = AlignIO.read(StringIO(">Alpha\nACGT\n>Gamma\nACGA"), "fasta")
two = AlignIO.read(StringIO(">Alpha\nAC--GT\n>Gamma\nAC--GA"), "fasta")
calc = DistanceCalculator("identity", skip_letters=("-", "*"))
print(calc.get_distance(one)["Alpha", "Gamma"])   # actual: 0.25
print(calc.get_distance(two)["Alpha", "Gamma"])   # actual: 0.5
# expected 0.25 in both cases: the same single substitution, and the extra
# columns are skipped gaps, not mismatches.

For reference, the same behaviour applies to Bio.Align.MultipleSeqAlignment and
Bio.Align.Alignment input, since get_distance() has a branch for each and both
call the same _pairwise().

A fix would make the identity branch accumulate max_score over the same
non-skipped positions that it uses for score, which also makes it consistent
with the scoring-matrix branch. The max_score == 0 guard further down already
handles the "everything was skipped" case, returning the maximum scaled distance
of 1.