Site WWW de Laurent Bloch
Slogan du site

ISSN 2271-3905
Cliquez ici si vous voulez visiter mon autre site, orienté vers des sujets moins techniques.

Pour recevoir (au plus une fois par semaine) les nouveautés de ce site, indiquez ici votre adresse électronique :

Un article des "Communications of the Association for Computing Machinery" (CACM)
La comparaison de séquences biologiques améliorée par l’algèbre
De la biologie pour les informaticiens
Article mis en ligne le 23 novembre 2016
dernière modification le 26 novembre 2016

par Laurent Bloch

Les numéros d’octobre et novembre 2016 des Communications of the Association for Computing Machinery (CACM) ont fait la part belle à la biologie, soit par les perspectives de renouvellement d’algorithmes utilisés quotidiennement par les biologistes, soit par une approche algorithmique de la théorie de l’évolution.

Séquences biologiques : résumé pour informaticien

L’exposé qui suit est très résumé, il n’a pour prétention que de donner une idée simplifiée des mécanismes étudiés par la biologie moléculaire.

Un des sports favoris des chercheurs en biologie moléculaire est la comparaison de séquences qui représentent des molécules biologiques, lesquelles viennent en trois variétés : protéines, acide désoxyribonucléique (ADN), acide ribonucléique (ARN).

Une protéine est une macromolécule formée d’une ou de plusieurs chaînes polypeptidiques. Chacune de ces chaînes est constituée de l’enchaînement d’acides aminés liés entre eux par des liaisons peptidiques. Les acides aminés, dans le monde vivant tel que nous le connaissons à la surface de la terre, sont au nombre de 20 (22 dans certains cas particuliers que nous omettons ici), et leur succession dans une chaîne polypeptidique est la séquence de la protéine (ou d’une sous-unité de protéine le cas échéant). Connaître la séquence d’une protéine ne livre pas toutes ses propriétés, mais beaucoup d’entre elles, notamment sa structure à trois dimensions dans l’espace, par calcul des repliements permis par l’agencement des acides aminés, ce qui fait l’objet d’étude de la modélisation moléculaire. On peut considérer les séquences protéiques comme des textes écrits dans un alphabet de 20 caractères, les acides aminés. Dans les bases de données de séquences protéiques chaque acide aminé est représenté conventionnellement par une lettre de l’alphabet.

Les macromolécules (on dit aussi polymères) d’ADN sont de longues chaînes de motifs moléculaires élémentaires, les nucléotides, dont il existe quatre variétés : l’adénine (notée A), la thymine (notée T), la cytosine (notée C) et la guanine (notée G). Le paradigme de la biologie moléculaire postule que l’information génétique est formulée par un texte, le génome, écrit dans un alphabet de quatre lettres, A, T, G, C, et que la connaissance de ce texte permet de connaître les fonctions de l’organisme considéré, sans avoir à entrer dans des considérations supplémentaires d’ordre physico-chimique.

L’ARN ressemble à l’ADN, mais à la place de la thymine (T) on y a l’uracile (U). L’ADN est très stable chimiquement, contrairement à l’ARN. De nombreux biologistes pensent que les premières manifestations de la vie sur terre reposaient sur l’ARN.

Pour qu’une protéine (alphabet de 20 caractères) soit synthétisée dans un organisme, il faut que la séquence correspondante existe dans son génome (alphabet de 4 caractères). Deux caractères (nucléotides) d’un texte d’ADN ne permettent que 16 combinaisons, pour traduire une séquence d’ADN en séquence protéique il faut donc trois nucléotides pour représenter un acide aminé, selon le code génétique. Trois nucléotides constituent un codon qui peut prendre 64 valeurs, de ce fait plusieurs configurations de nucléotides peuvent représenter le même acide aminé : on dit que le code génétique est dégénéré.

Comment cela se passe-t-il ? Des séquences codantes [1] de l’ADN génomique (les exons) sont transcrites en ARN dit messager par une enzyme, l’ARN-polymérase, nucléotide pour nucléotide, en remplaçant la thymine par l’uracile. Le brin d’ARN ainsi produit est ensuite traité par le ribosome, qui est une machinerie très complexe de protéines et d’ARN, qui va traduire le brin d’ARN messager en chaîne polypeptidique, pour finalement obtenir une protéine.

La transcription dépend des circonstances : notre organisme ne synthétise pas les mêmes protéines le matin que le soir, ni lorsque nous sommes en bonne santé que lorsque nous sommes malade, et bien sûr cela diffère d’un organe à l’autre de notre corps. Ce sont des gènes différents qui sont exprimés. L’étude de ces phénomènes permet de caractériser le transcriptome. Le génome humain comporte 3,4 milliards de nucléotides et quelques 25 000 gènes.

Le bon déroulement du processus suppose que la transcription et la traduction prennent bien les codons correctement, trois nucléotides par trois nucléotides, et s’arrêtent à un codon STOP (UAG, UGA ou UAA), c’est le cadre de lecture. Si, pour une raison ou une autre, le cadre de lecture se décale, par exemple parce qu’un nucléotide a été détruit ou au contraire inséré, le processus de traduction va donner des acides aminés complètement différents de ce que promettait le génome, le résultat se nomme frameshift mutation. La traduction peut aussi s’interrompre prématurément par l’occurrence erronée d’un codon STOP. De telles mutations peuvent être à l’origine de maladies génétiques telles que la drépanocytose et la thalassémie.

Analyse de séquences

La signification de cette locution est imprécise, mais lorsque j’ai débarqué sur la planète des biologistes en 1991 ils n’avaient qu’elle à la bouche. Ils commençaient à se résigner, non sans réticence, à l’idée que la recherche se fondait désormais sur l’analyse d’un texte. Il reste bien sûr d’autres champs de recherche en biologie, mais le champ de loin le plus productif est l’analyse du génome. Tous les trois ans je reçois la visite d’un de mes anciens étudiants du cours Pasteur d’informatique en biologie qui travaille maintenant principalement pour les National Institutes of Health (NIH) à Bethesda dans la banlieue de Washington : il y a des core labs dont la mission consiste à produire la séquence, et ensuite tout se passe exclusivement sur ordinateur. En France, dans ce domaine comme dans d’autres, les forces conservatrices ont encore beaucoup de pouvoir et elles réussissent à imposer des solutions bancales, en installant des chercheurs dans les laboratoires de séquençage, ce qui donne de la mauvaise recherche et du mauvais séquençage. Effectuer un séquençage de qualité n’est pas facile, et il est bien sûr indispensable que les ingénieurs qui s’en chargent soient en contact permanent avec la recherche en aval, mais tout mélanger ne donne rien de bon.

Dès qu’un biologiste a une séquence neuve, il n’a rien de plus pressé que de la comparer à d’autres séquences : séquences d’autres organismes de la même espèce, ou d’espèces différentes, pour repérer des similitudes ou des différences significatives. Cette comparaison devra tenir compte du fait que le monde du vivant comporte de l’imprévu : comme signalé plus haut, la transcription et la traduction de séquences peut comporter des erreurs, il y a aussi des polymorphismes, qui sont des variations locales du génome d’un individu à un autre. On ne cherchera donc pas une identité entre séquences ou parties de séquences, mais plutôt une plus ou moins grande ressemblance ; bref, c’est un problème d’optimisation : on attribue un coût aux différences (insertion, suppression, mutation) et on cherche parmi toutes les séquences à comparer celle dont l’alignement avec la nôtre donne le coût le plus faible.

Un autre article de ce site décrit l’algorithme de Needleman et Wunsch, typique de ce genre d’approche. L’idée est de « faire coulisser » les séquences l’une le long de l’autre pour trouver la position où la coïncidence est la meilleure, au besoin en insérant des espaces vides ici ou là pour que cela « colle mieux » (évidemment, tout est dans la bonne façon de faire coulisser et de repérer ce qui colle mieux, et d’en donner une définition formelle).

Ce type de démarche, si on l’effectue de manière naïve par force brute, conduit à effectuer de très nombreuses fois la même comparaison, ce qui donne une complexité exponentielle. La méthode classique pour échapper à cette malédiction est la programmation dynamique, qui consiste à enregistrer (ici dans un tableau) les résultats des comparaisons déjà effectuées pour ne plus avoir à les recommencer. Si les séquences à aligner sont de longueurs respectives n et m la programmation dynamique permet de passer d’une complexité en O(kn) à O(n.m), ce qui est remarquablement satisfaisant.

Supposons maintenant que nous disposions d’un ordinateur multiprocesseur, ou même simplement multicœur : chaque alignement de deux séquences est indépendant des autres alignements, il suffira, après les avoir tous calculés, de choisir le ou les meilleurs. Il suffit donc de partager le stock de séquences à comparer à la nôtre (nouvelle) en autant de tas qu’il y a de processeurs ou de cœurs dans notre système, et de répartir ces tas entre les processeurs disponibles, ce qui divisera le temps de calcul par le nombre de processeurs (enfin pas tout à fait, parce que le système d’exploitation consommera un peu de temps pour distribuer le travail aux processeurs).

Bon, la programmation dynamique est une idée judicieuse connue depuis des décennies, la répartition de calculs indépendants les uns des autres entre des processeurs distincts est une idée de bon sens, mais nous allons voir qu’il est possible de faire mieux, avec une approche plus subtile.

Parallélisation d’algorithmes de programmation dynamique

Si nous nous remémorons l’algorithme de Needleman et Wunsch nous voyons qu’il exploite le fait que pour calculer la valeur de la case Mi,j du tableau il faut (et il suffit de) connaître Mi-1,j, Mi,j-1 et Mi-1,j-1 et de prendre la plus grande de ces trois valeurs, ce que l’on peut représenter par un schéma ainsi :

C’est un grand progrès par rapport à l’algorithme par force brute, de complexité exponentielle. Nous pouvons bien sûr répartir les alignements de séquences sur les processeurs dont nous disposons, mais nous devons calculer les valeurs des cases de chaque tableau dans l’ordre des indices, séquentiellement.

Pour aller plus loin, les mathématiques tropicales

Ne serait-il pas possible de faire mieux ? Saeed Maleki, Madanlal Musuvathi et Todd Mytkowicz, de Microsoft Research à Redmond (WA), proposent dans un article intitulé Efficient Parallelization Using Rank Convergence in Dynamic Programming Algorithms du numéro d’octobre 2016 des CACM (pp. 85-92) une méthode générale de parallélisation à grain fin applicable à toute une classe d’algorithmes de programmation dynamique élégamment nommés linear-tropical dynamic programming problems (LTDP pour les intimes) [2] ; James Larus de l’École polytechnique fédérale de Lausanne en explique la portée dans son article introductif (p. 84), The Power of Parallelizing Computations.

Contrairement à ce qui est souvent le cas pour ce genre d’algorithme, la méthode remarquable de Messieurs Maleki, Musuvathi et Mytkowicz ne repose pas sur une habileté diabolique pour éliminer les cas redondants ou pour agencer les données, mais bien sur une démarche d’algèbre fondamentale. Au lieu de recourir à une parallélisation en front d’onde, qui calculerait tous les sous-problèmes d’une étape simultanément, et sans d’ailleurs que ce soit incompatible avec ce type d’optimisation, ils proposent de briser la dépendance entre les étapes successives en laissant dans le tableau des valeurs incorrectes, quitte à les rectifier ultérieurement.

Je serais bien incapable de discuter les détails de la méthode. Pour citer les auteurs : « On peut voir un calcul LTPD comme une séquence de multiplications matricielles dans un demi-anneau tropical où le demi-anneau est formé en prenant l’addition comme opération multiplicative et le maximum comme opération additive. Cet article démontre que plusieurs problèmes d’optimisation importants tels que Viterbi, LCS [longest common subsequence, NdT], Smith-Waterman [un proche cousin de Needleman-Wunsch, NdT] et Needleman-Wunsch (les deux derniers utilisés en bioinformatique pour l’alignement de séquences) sont LTPD. » L’article est assez compréhensible, j’y renvoie le lecteur pour la partie mathématique.

Mesures de performances, HPC vs systèmes distribués

La dernière section de l’article est consacrée à des mesures de performances. Les essais se sont déroulés sur des machines à mémoire partagée de type super-ordinateur (calcul haute performance, HPC) et sur des systèmes distribués, mais pour avoir les résultats respectifs sur ces deux types de systèmes il faudra se reporter au texte plus complet des Proceedings of the 19th ACM SIGPLAN symposium on Principles and practice of parallel programming, Parallelizing dynamic programming through rank convergence, où ce travail a été présenté initialement. Reste dans les CACM la comparaison entre les algorithmes de programmation dynamique classiques et celui de nos auteurs sur le système distribué Stampede du Texas Advanced Computing Center pour quatre problèmes LTPD : Viterbi, Smith-Waterman, Needleman-Wunsch et LCS, ce avec différentes tailles de paquets sur le réseau (1024, 2048, 4096, 8192 et 16 384) et différents nombres de processeurs (1 à 128 processeurs).

Les résultats varient selon la taille des paquets et le nombre de processeurs, mais dans la plupart des cas la version LTPD de l’algorithme procure une accélération d’un facteur compris entre 20 et 100.

Conclusion

Pour des problèmes que l’on pensait résolus depuis des années et pour lesquels n’étaient plus attendus que des progrès à la marge, une méthode entièrement fondée sur une approche très théorique obtient des améliorations très significatives. De surcroît cette méthode donne des résultats prometteurs avec des systèmes distribués à base d’ordinateurs (presque) ordinaires, ce qui réduit d’autant la nécessité du recours à des super-ordinateurs spécialisés et hors de prix, et ce pour des types de calculs effectués quotidiennement par de nombreuses équipes de recherche ainsi que dans l’industrie.