Distances évolutives

Prérequis

Cette section suppose acquise la construction d'un alignement multiple et la lecture d'un fichier FASTA, vues en 5.1 et en 1.3. Les matrices de substitution et la distinction entre score et distance sont introduites en 2.5. Il est utile de savoir qu'un alignement n'est pas une vérité biologique : les colonnes sont une hypothèse d'homologie, et les gaps y occupent une place particulière. Nous travaillerons ici sur de l'ADN, avec des séquences alignées de même longueur.

Introduction

Deux séquences homologues diffèrent à un certain nombre de positions. Compter ces différences paraît être l'opération la plus simple de la phylogénie moléculaire : sur un alignement de $L$ colonnes, on divise le nombre de colonnes discordantes par $L$ et l'on obtient un nombre entre 0 et 1. Ce nombre, appelé distance observée ou $p$-distance, a une vertu rare : il se calcule sans modèle et se vérifie à la main. Il a aussi un défaut majeur. Deux séquences peuvent avoir subi plusieurs substitutions au même site depuis leur ancêtre commun ; une partie de ces événements est invisible, et la différence observée sous-estime le nombre réel de substitutions.

L'enjeu biologique est considérable. Une distance mal corrigée fausse les branches d'un arbre, rapproche artificiellement des taxons à évolution rapide et peut faire croire à une divergence récente là où l'alignement ne contient plus de signal. Les virus, les bactéries et les gènes du développement illustrent toutes les situations où les taux de substitution varient fortement entre lignées. La composition en bases, le biais de transition et l'hétérogénéité des taux imposent de choisir un modèle explicite plutôt que d'appliquer une formule unique à tous les jeux de données.

L'enjeu computationnel est double. D'abord, il faut calculer une matrice complète en un temps raisonnable : $O(n^2)$ paires pour $n$ séquences, chacune en $O(L)$, soit $O(n^2L)$ opérations. Ensuite, il faut décider ce que l'on fait des sites ambigus, des gaps et des paires saturées. Une matrice qui contient des valeurs infinies, négatives ou silencieusement tronquées n'est pas une matrice de distances utilisable par Neighbor-Joining. Cette section construit donc des fonctions explicites, vérifiables et testées, avant de les appliquer à un alignement synthétique de démonstration.

Théorie

Distance observée et distance évolutive

Soient deux séquences alignées $s$ et $t$ de longueur $L$. On appelle distance observée, ou $p$-distance, la proportion de sites où les deux caractères diffèrent :

Dans cette section, $L$ désigne le nombre de sites comparables, c'est-à-dire les colonnes où les deux caractères sont A, C, G ou T. Les gaps et les codes d'ambiguïté sont exclus du dénominateur. Ce choix doit être documenté, car il change la valeur de $p$ : une colonne N contre A n'apporte aucune information sur la différence, mais elle ne doit pas être comptée comme une identité. Une séquence de 1000 pb avec 100 différences donne $p = 0{,}1$, soit 10 % de sites discordants.

La distance évolutive $d$ est une autre quantité. Elle représente le nombre attendu de substitutions par site depuis l'ancêtre commun, en comptant les substitutions multiples et les retours en arrière. Ce n'est pas une proportion observée : c'est un paramètre estimé sous un modèle. La relation entre $p$ et $d$ dépend des hypothèses sur les taux de substitution. Dire que $d = 0{,}1$ signifie qu'en moyenne 0,1 substitution par site s'est produite sur les deux branches, et non que 10 % des sites ont définitivement changé.

Substitutions multiples et modèle de Jukes-Cantor

Le modèle de Jukes et Cantor (1969) est le plus simple des modèles d'évolution moléculaire. Il suppose que les quatre bases ont la même fréquence (1/4 chacune), que toutes les substitutions sont également probables et que les sites évoluent indépendamment. Sous ces hypothèses, la probabilité qu'un site soit identique entre deux séquences séparées par une distance $d$ vaut :

La probabilité qu'il diffère est donc $p(d) = 1 - P(\text{identique})$, soit :

Cette formule montre le phénomène central de la section : $p$ n'est pas proportionnel à $d$. Pour les petites distances, $e^{-4d/3} \approx 1 - 4d/3$, donc $p \approx d$ : les substitutions multiples sont rares et la correction est faible. Pour les grandes distances, la courbe s'aplatit et $p$ tend vers $3/4$. Même si les séquences continuaient à diverger, la proportion de différences observées ne dépasserait pas cette limite en moyenne : les sites ont tellement changé qu'ils ne se souviennent plus de leur état initial.

En inversant la relation, on obtient la correction de Jukes-Cantor :

Le domaine de validité est $0 \leq p < 0{,}75$. À $p = 0{,}75$, l'argument du logarithme s'annule et $d$ diverge. Pour $p > 0{,}75$, la formule n'a pas de solution réelle : le modèle est saturé, ou ses hypothèses sont contredites par les données. Il ne faut jamais remplacer la valeur par un plafond arbitraire ; il faut signaler la paire et examiner l'alignement.

Exemple travaillé à la main : $p = 0{,}1$

Prenons un alignement de 1000 pb dans lequel 100 colonnes sont discordantes, soit $p = 0{,}1$. La distance observée vaut 0,1. Pour corriger, on calcule d'abord le facteur :

Le logarithme népérien de 0,866667 vaut environ $-0{,}143101$. En multipliant par $-3/4$ :

La distance corrigée est donc d'environ 0,1073 substitution par site. Sur 1000 pb, cela correspond à environ 107 substitutions attendues, contre 100 différences observées. La correction ajoute ici 7,3 % à la valeur brute. Elle est modeste parce que $p$ est petit ; pour $p = 0{,}5$, elle donne $d \approx 0{,}824$, soit une correction de 65 %. C'est précisément dans les branches longues que le choix du modèle devient décisif.

Vérifions la cohérence avec Kimura. Si les 100 différences se répartissent en 33 transitions et 67 transversions, on a $P = 0{,}033$ et $Q = 0{,}067$. La formule de Kimura donne alors $d \approx 0{,}1073$, presque la valeur de Jukes-Cantor : le rapport $P/Q = 1/2$ correspond au régime sans biais de transition. Nous revenons sur ce point ci-dessous.

Transitions, transversions et modèle de Kimura

Parmi les changements de base, on distingue deux classes. Une transition remplace une purine par une purine (A ↔ G) ou une pyrimidine par une pyrimidine (C ↔ T). Une transversion remplace une purine par une pyrimidine, ou l'inverse : A ↔ C, A ↔ T, G ↔ C, G ↔ T. Il y a donc deux types de transitions et quatre types de transversions, mais huit changements dirigés de transversion lorsque l'on tient compte du sens.

Le modèle de Kimura à deux paramètres (K2P, 1980) autorise des taux différents pour ces deux classes. Soient $P$ la proportion observée de sites portant une transition et $Q$ la proportion de sites portant une transversion, toujours rapportées au nombre de sites comparables. La distance corrigée s'écrit :

Les conditions de validité sont $1 - 2P - Q > 0$ et $1 - 2Q > 0$. Elles garantissent que les deux logarithmes sont définis. À la saturation, $P$ tend vers $1/4$ et $Q$ vers $1/2$ : les deux arguments tendent vers zéro et $d$ diverge. Le modèle ajoute donc un paramètre, le biais de transition, souvent résumé par le rapport $\kappa$ des taux de transition et de transversion. Ce rapport est estimé par maximisation de la vraisemblance dans un cadre phylogénétique ; il ne se lit pas directement comme le rapport $P/Q$ observé.

Le lien avec Jukes-Cantor est instructif. Lorsque $P/Q = 1/2$, c'est-à-dire lorsque les transitions représentent un tiers des différences, la formule de Kimura se réduit à celle de Jukes-Cantor. C'est le cas dans l'exemple précédent. À l'inverse, un fort biais de transition, fréquent chez les mammifères, produit des $P$ élevés et des branches plus longues que ne le laisserait croire la seule $p$-distance. Le choix entre JC et K2P n'est donc pas cosmétique : il dépend de la biologie du jeu de données.

GTR et la famille des modèles

Le modèle général à temps réversible (general time reversible, GTR) est le plus riche de la famille utilisée couramment pour l'ADN. Il possède six taux d'échange entre paires de bases, notés $r_{AC}, r_{AG}, r_{AT}, r_{CG}, r_{CT}, r_{GT}$, et quatre fréquences stationnaires $\pi_A, \pi_C, \pi_G, \pi_T$. Sa matrice de taux s'écrit, pour $i \neq j$ :

et la diagonale est choisie pour que chaque ligne somme à zéro. La réversibilité impose $\pi_i Q_{ij} = \pi_j Q_{ji}$ : le flux dans les deux sens est équilibré à l'équilibre. Après avoir fixé l'échelle globale des taux, les six échangeabilités ne comptent que pour cinq paramètres indépendants ; avec les trois fréquences indépendantes, GTR possède donc huit paramètres libres.

Les modèles plus simples sont des cas particuliers de GTR. Jukes-Cantor impose des fréquences égales et une seule échangeabilité. Kimura 2P conserve des fréquences égales mais distingue transitions et transversions. HKY ajoute le biais de fréquence à K2P. GTR généralise l'ensemble. Cette imbrication permet de comparer les modèles par des critères d'information (AIC, BIC) ou par des tests de rapport de vraisemblance, à condition de ne pas oublier la taille de l'échantillon. Un modèle plus complexe n'est pas automatiquement meilleur : avec un alignement court, il peut estimer des paramètres instables et produire des distances artificiellement précises.

Sous GTR, il n'existe pas de formule fermée simple pour la distance entre deux séquences. On estime $d$ par maximisation de la vraisemblance, numériquement, pour chaque paire. C'est plus coûteux que Jukes-Cantor, mais le coût reste modeste pour quelques dizaines de séquences. Pour un arbre complet, les programmes phylogénétiques utilisent souvent une matrice de distances pairwise comme point de départ, puis raffinent les paramètres globalement.

Saturation et perte de signal

La saturation est l'état dans lequel les substitutions multiples ont effacé une partie de l'information phylogénétique. Sur un graphique de $p$ en fonction de $d$, la courbe de Jukes-Cantor s'aplatit vers $3/4$ ; les transitions, souvent plus fréquentes, saturent avant les transversions. Avec des taux hétérogènes ou une composition biaisée, deux paires très divergentes peuvent donc avoir la même $p$-distance que des paires modérément divergentes, et même une $p$-distance supérieure à 0,75 si les séquences ne sont pas homologues.

Une saturation élevée n'est pas une simple curiosité théorique. Elle réduit la résolution des branches profondes, favorise les erreurs d'attraction des branches longues et rend les distances très sensibles au modèle. Les remèdes sont limités : allonger l'alignement, choisir un modèle mieux adapté, retirer les sites mal alignés, ou passer à une approche par vraisemblance qui modélise directement les colonnes. On ne corrige pas la saturation en multipliant une $p$-distance par un facteur arbitraire.

Distances entre séquences et entre profils

Jusqu'ici, chaque colonne de l'alignement était traitée comme un site net. Dans un alignement multiple, une colonne peut contenir des gaps, des ambiguïtés et des fréquences de bases. On parle alors de profil (profile) de colonne : un vecteur $f_X(a)$ qui donne la fréquence de chaque base $a$ dans la colonne $X$. Si l'on tire une base au hasard dans chaque profil, la probabilité qu'elles diffèrent est :

Cette quantité généralise la $p$-distance : lorsque chaque colonne ne contient qu'une seule base, $p_{XY}$ vaut 0 si les bases sont identiques et 1 sinon. On peut ensuite appliquer une correction de type Jukes-Cantor, mais l'hypothèse d'équilibre des fréquences doit être vérifiée. Pour comparer des profils plus riches, on utilise aussi des distances euclidiennes ou de Hellinger entre vecteurs de fréquences, ou des scores log-odds entre modèles de Markov cachés de profil. Ces scores ne sont pas toujours des distances métriques : il faut vérifier la symétrie et l'inégalité triangulaire avant de les donner à un algorithme qui les suppose.

Horloge moléculaire et limites

Une distance $d$ s'exprime en substitutions par site, pas en années. L'hypothèse d'horloge moléculaire stricte suppose un taux $\lambda$ constant par site et par unité de temps. Deux lignées séparées depuis $t$ unités de temps ont alors une distance attendue $d = 2\lambda t$, puisque chacune accumule $\lambda t$ substitutions. On peut donc estimer $t = d/(2\lambda)$ si $\lambda$ est calibré par un fossile, une date de divergence connue ou une série temporelle d'échantillons.

Cette hypothèse est souvent fausse. Les taux varient entre gènes, entre lignées, entre régions d'un génome et selon le mode de vie. Les modèles d'horloge relâchée autorisent ces variations, au prix de paramètres supplémentaires. Avant d'interpréter une distance comme un temps, il faut donc tester la constance des taux, documenter la calibration et distinguer un arbre de gènes d'un arbre d'espèces. Un arbre de distances peut être correct au sens statistique du modèle et faux biologiquement si l'échantillonnage, l'alignement ou l'homologie des sites sont inadéquats.

Implémentation

Les fonctions ci-dessous suivent une règle simple : toute valeur hors du domaine d'un modèle lève une exception explicite. On préfère un échec visible à une distance négative ou infinie qui contaminerait la matrice. Le code ignore les sites non ACGT dans $p$, conformément à la convention de la section, et conserve les noms des séquences dans l'ordre d'entrée.

Aucune sortie. Exécutez le code pour voir le résultat.

La sortie vérifiée le 2026-11-01 est :

0.1073
0.108
taxon_A [0.0, 0.0653, 0.0653, 0.1367]
taxon_B [0.0653, 0.0, 0.1367, 0.2158]
taxon_C [0.0653, 0.1367, 0.0, 0.2158]
taxon_D [0.1367, 0.2158, 0.2158, 0.0]

Trois points critiques méritent d'être soulignés. D'abord, la fonction p_distance normalise la casse mais ne supprime pas les gaps : elle les ignore site par site. Si l'on supprimait les colonnes contenant un gap avant le calcul, on changerait le dénominateur et donc la distance. Ensuite, distance_matrix ne remplit que le triangle supérieur puis recopie la valeur : la symétrie est garantie par construction, ce qui évite les matrices incohérentes. Enfin, la complexité est de $O(n^2L)$ en temps et $O(n^2)$ en mémoire. Pour quelques milliers de séquences, une implémentation vectorisée avec NumPy devient préférable, mais la logique reste la même.

La bibliothèque Biopython lit l'alignement et construit des arbres, mais elle ne fournit pas directement toutes les corrections de distance. On combine donc Bio.AlignIO pour l'entrée et nos fonctions pour les modèles. Cette séparation est utile : elle rend le calcul vérifiable et permet de comparer plusieurs modèles sur le même alignement sans changer le reste du pipeline.

Étude de cas

Nous utilisons l'alignement synthétique datasets/phylogeny/demo.fasta, créé pour ce cours. Il contient quatre séquences de 16 pb, nommées taxon_A à taxon_D. Il ne représente aucune espèce réelle et ne doit pas servir à une conclusion biologique. Sa petite taille permet en revanche de vérifier chaque distance à la main et d'observer l'effet de la correction.

La question est la suivante : quelles paires sont les plus proches, et la correction de Jukes-Cantor modifie-t-elle l'ordre des proximités ? Nous procédons en quatre étapes : lire l'alignement avec Biopython, vérifier les longueurs, calculer la matrice observée, puis la matrice corrigée.

Le script suivant est autonome et a été exécuté sur le fichier local. Il utilise AlignIO pour l'entrée et deux fonctions locales pour les distances.

import math
 
from Bio import AlignIO
 
alignment = AlignIO.read("datasets/phylogeny/demo.fasta", "fasta")
sequences = {record.id: str(record.seq).upper() for record in alignment}
 
 
def observed(seq1, seq2):
    return sum(left != right for left, right in zip(seq1, seq2)) / len(seq1)
 
 
def jc(p):
    return -0.75 * math.log(1.0 - 4.0 * p / 3.0)
 
 
names = list(sequences)
print("longueurs :", {name: len(seq) for name, seq in sequences.items()})
for title, transform in (("p-distance", None), ("Jukes-Cantor", jc)):
    print(title)
    print(" " * 10 + "".join(f"{name:>10}" for name in names))
    for name in names:
        values = [observed(sequences[name], sequences[other]) for other in names]
        if transform is not None:
            values = [0.0 if value == 0.0 else transform(value) for value in values]
        print(f"{name:>10}" + "".join(f"{value:10.4f}" for value in values))

La sortie vérifiée le 2026-11-01 est :

longueurs : {'taxon_A': 16, 'taxon_B': 16, 'taxon_C': 16, 'taxon_D': 16}
p-distance
             taxon_A   taxon_B   taxon_C   taxon_D
   taxon_A    0.0000    0.0625    0.0625    0.1250
   taxon_B    0.0625    0.0000    0.1250    0.1875
   taxon_C    0.0625    0.1250    0.0000    0.1875
   taxon_D    0.1250    0.1875    0.1875    0.0000
Jukes-Cantor
             taxon_A   taxon_B   taxon_C   taxon_D
   taxon_A    0.0000    0.0653    0.0653    0.1367
   taxon_B    0.0653    0.0000    0.1367    0.2158
   taxon_C    0.0653    0.1367    0.0000    0.2158
   taxon_D    0.1367    0.2158    0.2158    0.0000

Les valeurs observées sont des fractions simples : 1/16 pour A–B et A–C, 2/16 pour A–D et B–C, et 3/16 pour B–D et C–D. La correction augmente toutes les distances, mais pas dans la même proportion. Pour 1/16, elle passe de 0,0625 à 0,0653, soit +4,4 %. Pour 3/16, elle passe de 0,1875 à 0,2158, soit +15,1 %. L'ordre des proximités reste inchangé : A est plus proche de B et de C que B ne l'est de C, et D est le taxon le plus éloigné des trois autres. Cette stabilité est rassurante, mais elle ne prouve pas que l'arbre soit correct : avec 16 sites, l'incertitude d'échantillonnage est énorme.

Le modèle de Kimura donne des valeurs légèrement différentes pour certaines paires. La paire A–B, qui ne contient qu'une transversion, passe à 0,0657 ; la paire B–C, qui contient deux transversions, passe à 0,1387. La paire B–D conserve 0,2158, car son rapport transitions/transversions correspond au régime de Jukes-Cantor. Aucune paire n'est saturée : la plus grande $p$-distance vaut 0,1875, très loin de 0,75. Le jeu est donc utilisable pour illustrer la construction d'un arbre en 5.3, mais il ne faut pas lui demander une phylogénie réelle.

Erreurs fréquentes

Quiz

Exercices

python -m pytest exercises/module-5/5.2-distances-evolutives -q

python -m pytest exercises/module-5/5.2-distances-evolutives -q

python -m pytest exercises/module-5/5.2-distances-evolutives/projet -q

Résumé

• La $p$-distance est la proportion observée de sites différents ; elle dépend de l'alignement et du traitement des gaps.
• Les substitutions multiples rendent $p$ inférieur à $d$ : la correction de Jukes-Cantor vaut $d = -3/4 \ln(1 - 4p/3)$.
• Le modèle de Kimura distingue transitions et transversions et s'écrit avec deux proportions $P$ et $Q$.
• GTR combine six échangeabilités et quatre fréquences ; JC, K2P et HKY sont des cas particuliers emboîtés.
• La saturation correspond à une perte de signal ; à $p = 0{,}75$, la correction de JC diverge.
• Une matrice de distances doit être symétrique, à diagonale nulle, et documenter les paires non calculables.
• Une distance est en substitutions par site, pas en années : l'horloge moléculaire est une hypothèse à tester.
• Le choix du modèle et le trimming changent les valeurs ; il faut les rapporter avec les résultats.

Références

• Jukes, T. H. & Cantor, C. R. (1969). Evolution of protein molecules. In Mammalian Protein Metabolism, vol. 3, p. 21–132. Academic Press.
• Kimura, M. (1980). A simple method for estimating evolutionary rates of base substitutions through comparative studies of nucleotide sequences. Journal of Molecular Evolution, 16(2), 111–120.
• Saitou, N. & Nei, M. (1987). The neighbor-joining method: a new method for reconstructing phylogenetic trees. Molecular Biology and Evolution, 4(4), 406–425.
• Felsenstein, J. (1985). Confidence limits on phylogenies: an approach using the bootstrap. Evolution, 39(4), 783–791.
• Capella-Gutiérrez, S., Silla-Martínez, J. M. & Gabaldón, T. (2009). trimAl: a tool for automated alignment trimming in large-scale phylogenetic analyses. Bioinformatics, 25(15), 1972–1973.
• Cock, P. J. A. et al. (2009). Biopython: freely available Python tools for computational molecular biology and bioinformatics. Bioinformatics, 25(11), 1422–1423.
• Documentation Biopython 1.88. Bio.Phylo.TreeConstruction, Bio.AlignIO. https://biopython.org/docs/1.88/api/Bio.Phylo.TreeConstruction.html