L'alignement multiple en étoile ⭐

Introduction

Dans le domaine de la bio­in­for­ma­tique l'ali­gne­ment de séquences est une opé­ra­tion cou­rante pour com­pa­rer des séquences bio­lo­giques (ADN, ARN, pro­téines). Nous avons vu dans un pré­cé­dent article le prin­cipe de l'algorithme de Need­le­man-Wunsch qui pro­duit un ali­gne­ment glo­bal opti­mal de deux séquences.

Au-delà de l'alignement par paire, l'ali­gne­ment mul­tiple de séquences (Mul­tiple Sequence Ali­gn­ment ou MSA) per­met de com­pa­rer simul­ta­né­ment plu­sieurs séquences afin de mettre en évi­dence des régions conser­vées, de repé­rer des motifs fonc­tion­nels ou de faci­li­ter des ana­lyses phy­lo­gé­né­tiques. La géné­ra­li­sa­tion naïve de l'algorithme de Need­le­man-Wunsch pour plus de deux séquences est mal­heu­reu­se­ment irréa­li­sable à cause de sa com­plexi­té tem­po­relle en (pour séquences de nucléo­tides).

Dans ce billet, je vais vous pré­sen­ter l'intérêt et le prin­cipe de l'ali­gne­ment mul­tiple en étoile (Star Ali­gn­ment ou par­fois Cen­ter Star Ali­gn­ment), une approche qui repose sur la réa­li­sa­tion d'alignements par paires, puis leur fusion autour d'une séquence cen­trale de réfé­rence.

1. Exemple d'incohérences entre alignements de paires

L'une des rai­sons qui rend le MSA si com­plexe est le fait que les ali­gne­ments de paires de séquences, bien qu'optimaux lorsqu'on les consi­dère indi­vi­duel­le­ment, sont sou­vent inco­hé­rents entre eux. Il n'est donc pas pos­sible de fusion­ner sim­ple­ment ces ali­gne­ments pour pro­duire un ali­gne­ment mul­tiple.

Pre­nons un exemple simple de 3 séquences de 10 nucléo­tides pour mon­trer com­ment se mani­festent ces inco­hé­rences. Par com­mo­di­té je choi­sis les para­mètres d'alignement sui­vants : +1 pour une iden­ti­té (match), -1 pour une sub­sti­tu­tion (mis­match) et -1 éga­le­ment pour une lacune (gap). Avec les 3 séquences (seq1, seq2, seq3) sui­vantes nous obte­nons 3 paires d'alignements (ali1/​2, ali1/​3, ali2/​3).

Pour faci­li­ter la com­pa­rai­son entre ces 3 ali­gne­ments et mettre en évi­dence leurs inco­hé­rences on peut numé­ro­ter les nucléo­tides. Ain­si :

• D'après l'alignement 1, le nucléo­tide 1 de la séquence 1 (le deuxième A) est homo­logue au nucléo­tide 2 de la séquence 2.
• D'après l'alignement 2, ce même nucléo­tide 1 de la séquence 1 est homo­logue au nucléo­tide 0 de la séquence 3.

On devrait donc en déduire que le nucléo­tide 2 de la séquence 2 et le nucléo­tide 0 de la séquence 3 sont homo­logues. L'homologie est en effet une rela­tion tran­si­tive. Mais nous consta­tons au contraire que dans l'alignement 3 le nucléo­tide 2 de la séquence 2 n'est homo­logue à rien du tout dans la séquence 3, et que le nucléo­tide 0 de la séquence 3 est homo­logue au nucléo­tide 0 de la séquence 2.

Les 3 ali­gne­ments pris ensemble sont donc inco­hé­rents, on ne peut pas les fusion­ner en un ali­gne­ment mul­tiple unique sans contre­dire cer­taines por­tions de ces ali­gne­ments par paires.

2. Principe de l'alignement en étoile

Dans l'alignement en étoile on com­mence par choi­sir l'une des séquences qu'on dési­gne­ra comme cen­trale et qui ser­vi­ra de réfé­rence. On aligne cette séquence suc­ces­si­ve­ment avec toutes les autres, mais on n'aligne pas ces autres séquences entre elles.

Sché­ma­ti­que­ment, on peut donc pla­cer la séquence cen­trale au milieu d'une étoile dont toutes les branches repré­sentent un ali­gne­ment par paire avec l'une des autres séquences (voir sché­ma ci-des­sous). Le fait de ne pas consi­dé­rer les ali­gne­ments entre les séquences de la péri­phé­rie de l'étoile per­met d'éliminer les inco­hé­rences, mais on perd aus­si bien sûr de l'information.

3. Choisir une séquence centrale

Il existe plu­sieurs stra­té­gies pour choi­sir la séquence cen­trale de manière à mini­mi­ser cette perte d'information. La plus com­mune est de choi­sir la séquence qui est en moyenne la plus proche de toutes les autres. On espère alors de cette manière pri­vi­lé­gier les ali­gne­ments par paire conte­nant le moins d'erreurs d'alignement.

Plu­sieurs détails impor­tants peuvent varier dans la mise en œuvre de cette stra­té­gie. Par exemple quelle for­mule on décide d'appliquer pour cal­cu­ler les dis­tances entre les séquences. Péna­lise-t-on ou ignore-t-on les gaps ? Est-ce qu'on uti­lise une dis­tance brute comme la p-dis­tance ou est-ce qu'on applique une cor­rec­tion comme avec la dis­tance de Jukes-Can­tor ? Com­ment cal­cule-t-on cette moyenne des dis­tances ? Avec ou sans pon­dé­ra­tion ? Ou peut-être pré­fè­re­ra-t-on uti­li­ser la médiane fina­le­ment ?

Bref, res­tons simple et pour­sui­vons avec notre exemple. L'alignement 1 montre que les séquences 1 et 2 ont une p-dis­tance nulle car il n'y a que des matchs et aucun mis­matchs. Le cal­cul de la p-dis­tance consiste en effet à divi­ser le nombre de nucléo­tides dif­fé­rents par le nombre total de nucléo­tides com­pa­rés, ce qui implique d'ignorer toutes les posi­tions avec des gaps (le "p" signi­fie "pro­por­tion").

L'alignement 2 montre un seul mis­match sur 7 posi­tions sans gap, donc la p-dis­tance est de 1/​7 entre les séquences 1 et 3. Enfin l'alignement 3 montre 3 mis­matchs sur 8 posi­tions sans gap, la dis­tance entre les séquences 2 et 3 est donc de 3/​8. Cela nous per­met de cal­cu­ler la dis­tance moyenne d'une séquence aux deux autres séquences :

La séquence 1 est celle qui appa­rait en moyenne la plus proche, c'est donc elle qu'on choi­sit ici comme séquence cen­trale. On rejette donc com­plè­te­ment l'alignement 3 qui n'inclut pas la séquence cen­trale.

4. Construire une séquence centrale cohérente

Dans une implé­men­ta­tion simple du prin­cipe de l'alignement en étoile, l'étape qui suit le choix de la séquence cen­trale est celui de l'agglomération de toutes les ver­sions de cette séquence pour construire une unique ver­sion de réfé­rence. Cette ver­sion finale doit incor­po­rer toutes les lacunes (gaps) des ali­gne­ments par paires rete­nus. Voi­ci les deux ver­sions de la séquence cen­trale de notre exemple pro­ve­nant res­pec­ti­ve­ment des ali­gne­ments 1 et 2 :

Pour agglo­mé­rer deux ver­sions il suf­fit de les par­cou­rir simul­ta­né­ment avec deux poin­teurs dif­fé­rents. Ici je parle de poin­teur par abus de lan­gage pour dési­gner sim­ple­ment une variable qui contient la posi­tion qu'on lit dans une chaine de carac­tère (on peut accé­lé­rer l'algorithme si on uti­lise un vrai poin­teur qui contient l'adresse mémoire de cette posi­tion mais je ne ren­tre­rai pas dans ces détails d'optimisation). On pro­cède alors comme suit :

• Lorsqu'on lit avec les deux poin­teurs un même carac­tère non-gap alors on ajoute ce carac­tère à la nou­velle ver­sion et on incré­mente les deux poin­teurs. Remar­quez que ce carac­tère non-gap est néces­sai­re­ment le même dans les deux séquences, si ce n'est pas le cas c'est qu'il y a une erreur dans l'implémentation de l'algorithme.
• Lorsqu'on lit un gap avec les deux poin­teurs alors on incré­mente de la même façon les deux poin­teurs et on ajoute un gap à la nou­velle ver­sion.
• Lorsqu'on lit avec un seul des deux poin­teurs un gap alors on ajoute ce gap à la nou­velle ver­sion et on incré­mente seule­ment ce poin­teur, mais pas l'autre.

Voi­ci ce que cela donne sur notre exemple :

• Le poin­teur p1 et le poin­teur p2 sont ini­tia­li­sés à 0 et lisent res­pec­ti­ve­ment "-" et "A". Seul le poin­teur p1 est donc incré­men­té et on ajoute un gap à la nou­velle ver­sion.
  
• Le poin­teur p1 vaut 1 et p2 vaut tou­jours 0, les deux lisent "A", donc "A" est ajou­té à la nou­velle ver­sion et les deux sont incré­men­tés.
  
• De la même manière on ajoute "A", "T", "A", "G". On a alors p1 = 6 et p2 = 5, p1 lit "-" et p2 lit "G", on ajoute alors un gap à la nou­velle ver­sion et seul p1 est incré­men­té à 7.
  
• Avec p1 = 7 et p2 = 5 on lit res­pec­ti­ve­ment encore "-" et "G", donc on ajoute à nou­veau un gap et on incré­mente seule­ment p1.
• Avec p1 = 8 et p2 = 5 on lit deux "G", donc on l'ajoute et on incré­mente les deux. Avec p1 = 9 et p2 = 6 on lit encore deux "G" donc on fait de même. Ensuite avec p1 = 10 et p2 = 7 on lit deux "C", donc on l'ajoute et on incré­mente les deux.
  
• Avec p1 = 11 et p2 = 8 on lit res­pec­ti­ve­ment "T" et "-", donc cette fois c'est seule­ment p2 qu'on incré­mente à 9 et on ajoute un gap.
  
• Je vous laisse finir ?

Voi­ci le résul­tat qu'on obtient, vous pou­vez consta­ter que cet algo­rithme per­met bien d'incorporer les gaps des deux ver­sions dans une seule nou­velle ver­sion de la séquence :

Si on a plus de deux ver­sions à fusion­ner, un algo­rithme naïf consiste à tout sim­ple­ment fusion­ner ces ver­sions pro­gres­si­ve­ment dans une boucle. L'ordre des fusions n'a aucune influence sur le résul­tat final.

5. Mettre en cohérence toutes les autres séquences

À cette étape nous connais­sons la ver­sion finale de la séquence cen­trale et nous avons aus­si dif­fé­rentes ver­sions locales de cette séquence cen­trale ali­gnées avec cha­cune des autres séquences. Ce sont ces autres séquences que nous vou­lons main­te­nant ali­gner avec la ver­sion finale de la séquence cen­trale.

Le prin­cipe est le même que dans l'étape pré­cé­dente. Les gaps man­quants dans ces autres séquences sont les mêmes que ceux qui manquent si on com­pare la ver­sion finale de la séquence cen­trale avec sa ver­sion locale pour une paire don­née. Il suf­fit donc de par­cou­rir chaque ali­gne­ment et de com­pa­rer les deux ver­sions de la séquence cen­trale : la finale et la locale.

• Lorsque les deux ver­sions de la séquence cen­trale ont le même carac­tère (gap ou non) alors on ne modi­fie pas la séquence à ali­gner, puis on incré­mente les deux poin­teurs de lec­ture.
• Lorsque les deux ver­sions sont dif­fé­rentes alors on ajoute un gap à cette posi­tion dans la séquence à ali­gner, puis on incré­mente seule­ment le poin­teur de lec­ture de la ver­sion finale de la séquence cen­trale.

Avec cette pro­cé­dure on obtient l'alignement mul­tiple final sui­vant :

Comme atten­du, cet ali­gne­ment mul­tiple est par­fai­te­ment cohé­rent avec les deux ali­gne­ments par paire ali1/​2 et ali1/​3, mais n'est que par­tiel­le­ment cohé­rent avec l'alignement ali2/​3. En fait si on com­pare ces ali­gne­ments plus pré­ci­sé­ment, on ne retrouve en com­mun dans l'alignement mul­tiple final et dans l'alignement ali2/​3 que le pre­mier "T" (nucléo­tide d'indice 3 de la séquence 2 et nucléo­tide d'indice 1 de la séquence 3), aucunes des autres asso­cia­tions ne sont conser­vées.

6. Évaluation de la qualité de l'alignement multiple

Comme on l'a vu, le fait d'avoir choi­si la séquence 1 comme séquence cen­trale nous a obli­gé à igno­rer les infor­ma­tions de l'alignement ali2/​3. Si on avait choi­si la séquence 2 ou bien la séquence 3 comme séquence cen­trale alors nous aurions igno­ré à la place d'autres infor­ma­tions et nous aurions obte­nu les ali­gne­ments mul­tiples sui­vants.

Ali­gne­ment mul­tiple en étoile avec la séquence 2 au centre :

Ali­gne­ment mul­tiple en étoile avec la séquence 3 au centre :

Com­ment savoir si un de ces ali­gne­ments ne serait pas meilleur que celui que nous avons obte­nu ?

Il existe plu­sieurs manières d'évaluer la qua­li­té d'un ali­gne­ment mul­tiple. On fait géné­ra­le­ment appel au cal­cul d'une fonc­tion objec­tive qu'on cherche à mini­mi­ser ou à maxi­mi­ser selon les cas. Je rap­pelle que le but est d'obtenir un maxi­mum de vraies posi­tions homo­logues et un mini­mum d'erreurs. Il est donc natu­rel d'extrapoler l'utilisation du score d'alignement qui a fait le suc­cès de l'algorithme de Need­le­man-Wunsch pour infé­rer cor­rec­te­ment les posi­tions homo­logues d'une paire de séquences.

On peut donc cal­cu­ler un score d'alignement mul­tiple en fai­sant la somme des scores d'alignement pour chaque paire de cet ali­gne­ment mul­tiple. Cette fonc­tion objec­tive s'appelle sum-of-pairs ou "somme des paires". On doit bien sûr igno­rer lors de ce cal­cul les paires de posi­tions avec deux gaps.

Com­men­çons par cal­cu­ler le score pour un ali­gne­ment simple :

Sur la ligne d'évaluation j'ai mis des points . pour repré­sen­ter les mis­matchs et les gaps (par sim­pli­ci­té tous les deux ont ici un score de -1), les plus + repré­sentent les matchs (ici avec un score de +1), et j'ai mis des zéros 0 pour deux gaps. On compte 7 plus et 6 points, donc le score de cet ali­gne­ment est de 7 - 6 = 1.

On fait la même chose pour les deux autres paires de l'alignement mul­tiple A (celui avec la séquence 1 cen­trale) on obtient de cette manière le score de -10 :

Sans sur­prise c'est la paire 2/​3 qui reçoit le score le plus bas puisque ses infor­ma­tions ont été igno­rées. Les ali­gne­ments B et C reçoivent res­pec­ti­ve­ment les scores -12 et -13 :

Ces scores plus bas nous montrent que ces ali­gne­ments alter­na­tifs ne sont pas meilleurs et que nous avons eu rai­son de choi­sir la séquence 1 comme séquence cen­trale.

Notez bien que des choix de para­mètres dif­fé­rents pour les matchs, les mis­matchs et les gaps chan­ge­raient les ali­gne­ments par paires et l'alignement mul­tiple final en modi­fiant les cal­culs des scores. On peut de plus choi­sir une autre fonc­tion objec­tive que la simple sum-of-pairs, à com­men­cer par la weigh­ted sum-of-pairs qui ajoute une pon­dé­ra­tion plus ou moins com­plexe dans le cal­cul de cette somme.

7. Intérêts et limites

L'alignement mul­tiple en étoile est simple à implé­men­ter et très rapide à mettre en œuvre. En fait ce sont les ali­gne­ments par paires qui consti­tuent le gou­lot d'étranglement prin­ci­pal de l'algorithme. Mais on peut évi­ter d'avoir à cal­cu­ler les ali­gne­ments de toutes les paires si notre stra­té­gie de choix de la séquence cen­trale ne le requiert pas (je ne détaille­rai pas ici les tech­niques de cal­culs de dis­tances entre séquences basées sur les k-mères qui ne néces­sitent pas d'alignements).

Mal­heu­reu­se­ment il n'existe pas de stra­té­gie de choix de la séquence cen­trale qui nous assure d'obtenir le meilleur ali­gne­ment mul­tiple. On pour­rait être ten­té d'essayer toutes les séquences cen­trales pos­sibles et de choi­sir l'alignement final ayant le meilleur score, mais cette méthode est beau­coup plus lente car elle nous oblige à cal­cu­ler tous les ali­gne­ments par paires.

Même avec une stra­té­gie aus­si cou­teuse l'alignement en étoile est peu per­for­mant avec des séquences assez éloi­gnées. Les résul­tats ne sont fiables que pour des séquences très proches. Cela est dû au fait que l'algorithme ignore une grande quan­ti­té d'informations dans les ali­gne­ments par paires qui ont été reje­tés. Ain­si, une erreur d'alignement com­mise dans l'un des ali­gne­ments par paire rete­nu ne peut pas être cor­ri­gée par les infor­ma­tions conte­nues dans les ali­gne­ments par paires reje­tés.

Des stra­té­gies d'alignement mul­tiple plus com­plexes et plus pré­cises ont donc été déve­lop­pées pour exploi­ter au maxi­mum l'ensemble des infor­ma­tions dis­po­nibles. Les ali­gne­ments mul­tiples basés sur la cohé­rence (consis­ten­cy based MSA) essaient par exemple de la cohé­rence des homo­lo­gies entre tous les ali­gne­ments par paires. Les ali­gne­ments mul­tiples pro­gres­sifs (pro­gres­sive MSA) essaient plu­tôt de construire un ali­gne­ment en ajou­tant les séquences une par une et en l'alignant de manière à prendre en compte les infor­ma­tions dans toutes les séquences déjà ali­gnées.

Références

• Gus­field, D. (1997). Mul­tiple string com­pa­ri­son - the Holy Grail. Algo­rithms on Strings, Trees and Sequences, 332-366.
• Gus­field, D. (1993). Effi­cient methods for mul­tiple sequence ali­gn­ment with gua­ran­teed error bounds. Bul­le­tin of mathe­ma­ti­cal bio­lo­gy, 55(1), 141-154.
• Gus­field, D. (1991). Effi­cient algo­rithms for infer­ring evo­lu­tio­na­ry trees. Net­works, 21(1), 19-28.

Annexe : implémentation de l'alignement en étoile en Python

L'algorithme d'alignement en étoile sui­vant n'est pas du tout opti­mi­sé. Cette manière d'écrire le code per­met de sépa­rer clai­re­ment les étapes pour que l'ensemble soit faci­le­ment com­pré­hen­sible, mais cela pro­voque des redon­dances. Réfé­rez-vous à l'article pré­cé­dent pour l'implémentation de l'algorithme de Needl­man-Wunsch qui pro­duit les ali­gne­ments par paires.

Mer­ci à Léo­pold Car­ron et Lee Mariault pour la relec­ture.