Détermination du génome coeur de 7 souches de la bactérie Escherichia coli

Nous avons vu en cours que la taille et le nombre de gènes entre génomes de différentes souches d'une même espèce bactérienne pouvaient varier. Le génome coeur (ou core genome) a été défini comme correspondant à l'ensemble des gènes communs à toutes les souches d'une même espèce (ensemble des gènes orthologues).

Comment pratiquement et rapidement peut-on identifier les gènes orthologues entre deux espèces ?

La technique la plus couramment utilisée consiste à analyser les résultats d'un programme appelé Blast qui permet de comparer rapidement une séquence donnée (sequence query) à un ensemble de séquences constituant une banque. Le résultat du Blast est un alignement local indiquant les régions présentant des similarités entre les deux séquences. Cette similarité est quantifiée par le pourcentage d'identité trouvé entre les deux séquences alignées, le score de l'alignement (qui a été normalisé et est appelé bit score) et la e-value qui correspond au nombre attendu d’alignement qui par chance aurait un score >= au score obtenu avec nos deux "vraies" séquences. Ce nombre  est évalué en regardant les scores des alignements générés par comparaison de séquences aléatoires ayant même longueur et même composition que la séquence requête. Plus la e-value est faible plus l'alignement obtenu est significatif.

Pour rechercher les gènes orthologues, la comparaison se fait souvent au niveau de leurs produits protéiques. Nous utiliserons donc ici, par abus de langage, les termes d'orthologues et de paralogues pour des séquences protéiques, ceci pour simplifier l'écriture.

Dans les fichiers joints vous trouverez donc le résultat de la comparaison de chaque protéine d'un génome (chaque séquence sera à son tour séquence query) avec l'ensemble des protéines d'un autre génome constituant donc la banque de séquences. Le programme utilisé pour la comparaison est BlastP (version du blast pour comparer une séquence protéique à une banque protéique).

Nomenclature utilisée et exemple
EcolB_EcolD.blt : chaque séquence du génome de la souche EcolB01 sera une séquence query qui sera comparée à l'ensemble des séquences protéiques de la souche EcolD01. Comme résultats on ne gardera que le premier hit (donc la séquence de EcolD01 qui possède la plus petite e-value). On parlera alors de Best Hit.
Si on prend l'exemple de la séquence EcolB01.ECS002 son best hit dans EcolD01 est la séquence EcolD01.THRA. On notera BH(EcolB01.ECS002, EcolD01.THRA) la e-value est de 0 et le score de 1590.
EcolD_EcolB.blt : chaque séquence du génome de la souche EcolD01 sera une séquence query qui sera comparée à l'ensemble des séquences protéiques de la souche EcolB0101. Comme résultats on ne gardera que le premier hit.
Si on prend l'exemple de la séquence EcolD01.THRA, son best hit est la séquence EcolB01.ECS002. On a donc BH(EcolD01.THRA,EcolB01.ECS002)

Dans cet exemple on a donc BH(EcolB01.ECS002, EcolD01.THRA) = BH(EcolD01.THRA,EcolB01.ECS002), les deux séquences sont réciproquement best hit l'une de l'autre. Dans ce cas on parle de Best Bidirectional Hits appelé dans la littérature BBH.

La première étape condition pour que deux protéines soient orthologues est qu'elles soient BBH. Ensuite, pour que cette relation soit un-un, il ne faut pas qu'une des deux séquences possède une séquence paralogue dont le score d'alignement serait supérieur à celui du BBH.

Nous devons donc vérifier qu'il n'y a pas eu de duplication du gène dans un des deux génomes qui serait récente et résulterait dans la présence d'une séquence paralogue (dans le même génome) ayant un score de d'alignement de Blast supérieur à celui obtenu entre les séquences BBH de nos deux génomes. Pour cela, pour chaque génome, chacune de ces séquences protéiques a été comparée à l'ensemble des autres séquences protéiques et le premier et donc meilleur hit obtenu est appelé Best Paralog, noté BP.
Les résultats de ces comparaisons sont trouvés dans les fichiers portant deux fois le nom de la souche.
Si nous reprenons l'exemple précédant :
EcolB_EcolB.blt : chaque séquence protéique de la souche EcolB01 est comparée à toutes les autres protéines de cette souche. Le best hit de la séquence EcolB01.ECS002 est EcolB01.ECS4869. Ces deux séquences sont donc best paralog et le score de l'alignement obtenu est 341.
EcolD_EcolD.blt : chaque séquence protéique de la souche EcolD01 est comparée à toutes les autres protéines de cette souche. Le best hit de la séquence EcolD01.THRA est EcolD01.METL. Ces deux séquences sont donc best paralog et le score de l'alignement obtenu est 344.


Les deux scores BP(EcolB01.ECS002,EcolB01.ECS4869) = 341 et BP(EcolD01.THRA,EcolD01.METL) = 344 < BBH((EcolB01.ECS002, EcolD01.THRA) = 1590

Il n'existe donc pas de séquences paralogues dans aucun des deux génomes possédant un score d'alignement supérieur au BBH. Donc nous pouvons conclure que EcolB01.ECS002 et EcolD01.THRA sont deux séquences orthologues (un-un).

Les données fournies :  (archive ci-jointe)

7 souches ont été analysées : EcolB, EcolD, EcolL,EcolM, EcolR, EcolU et EcolX.

Pour chacune d'entre elle il y a 7 fichiers : Comparaison de chaque protéine de la souche avec les autres protéines de la même souche (fichier ou le nom de la souche apparaît deux fois), comparaison de chaque protéine de la souche avec l'ensemble des protéines des 6 autres souches (dans le nom du fichier, le nom de la première souche est celle dont les séquences sont utilisées comme query et le deuxième est le nom de la souche dont les séquences ont servi à constituer la banque).
Dans chaque cas seul le best hit est reporté.
Le fichier contient le résultat du BlastP sous forme de tableau avec les colonnes correspondant aux données suivantes:

nom séquence query , nom séquence banque, % d'identité, longueur de l'alignement, nombre de substitutions, nombre de gaps (indels), position début alignement séquence query, position fin alignement séquence query, position début alignement séquence banque, position fin alignement séquence banque, e-value, score normalisé (bit score).

Pour trouver le génome coeur il va donc falloir comparer ces 7 fichiers et trouver les protéines orthologues dans les 7 génomes. Dans l'étape précédente, on a trouvé que les deux séquences EcolB01.ECS002 et EcolD01.THRA étaient orthologues. Maintenant, il va falloir rechercher le BBH de chacune d'entre-elles dans les autres génomes, vérifier que la séquence identifiée comme BBH dans le nouveau génome que l'on parte de EcolB01.ECS002 ou de EcolD01.THRA est bien la même et qu'il n'y a pas de paralogues dans ce nouveau génome de plus fort score blast.
Exemple si on regarde la souche EcolL.
On a BBH(EcolB01.ECS002,EcolL01.THRA) = 1595 et on a BBH(EcolD01.THRA,EcolL01.THRA). Le best paralog est BP(EcolL01.THRA, EcolL01.METL) = 344. Donc dans les trois souches pas de paralogues de score plus élevé. Donc ici les 3 protéines sont orthologues.

On fait de même avec les quatre génomes restants. Si on a trouvé une protéine orthologue dans chacun des génomes, alors la protéine appartient au génome coeur.

Le projet demandé
Faire le programme qui permet d'extraire le génome coeur en appliquant les principes énoncés ci-dessus et donner le résultat obtenu sur les 7 génomes analysés, c'est à dire le nombre et la liste des protéines constituant ce génome coeur.