Compareads et Commet

Après vous avoir par­lé de Méta­gé­no­mique et de filtre de Bloom, voi­ci un article fusion­nant les deux concepts en un algo­rithme !

Pour rap­pel, la méta­gé­no­mique vise à étu­dier le conte­nu géné­tique et géno­mique d'un milieu don­né. On passe géné­ra­le­ment par un séquen­çage NGS et on récu­père un grand nombre de courtes séquences d'ADN (lec­tures, ou "reads") pour chaque échan­tillon bio­lo­gique ana­ly­sé.

Tou­jours pour rap­pel, le filtre de Bloom est une struc­ture de don­nées sim­pliste qui sert à tes­ter si un élé­ment fait par­tie d'un ensemble de don­nées ou non. Cette struc­ture a l'avantage d'être très éco­nome en mémoire et très rapide en exé­cu­tion. Par contre, il est pos­sible d'avoir des faux-posi­tifs.

Alors, com­ment fusion­ner ces deux concepts ?

On peut ima­gi­ner indexer toutes les séquences d'un échan­tillon méta­gé­no­mique A dans le filtre de Bloom, puis tes­ter chaque séquence d'un échan­tillon B pour savoir si elle est pré­sente dans A ou non ! On obtient ain­si la liste et le nombre de séquences simi­laires entre deux échan­tillons méta­gé­no­mique et on parle alors de méta­gé­no­mique com­pa­ra­tive.

Plus exac­te­ment on ne va pas indexer les séquences, mais indexer tous les k-mers de chaque séquence. Un k-mer est une sous séquence de taille k. Par exemple la séquence GTAATAC se décom­pose en six 3‑mers : GTA, TAA, AAT, ATA, TAC. Elle se décom­pose aus­si en trois 5‑mers : GTAAT, TAATA et AATAC. Et ain­si de suite pour des 1‑mers, 2‑mers etc.

On n'indexe pas direc­te­ment les séquences mais les k-mers pour pou­voir iden­ti­fier des séquences "proches" et pas com­plè­te­ment iden­tiques. Si j'indexe GTAATAC dans le filtre de Bloom, et qu'ensuite je teste si TAATACT y est, la réponse sera non ; avec des séquences de 100 paires de bases, il n'y aura donc aucune séquence com­munes ente les deux jeux de don­nées.

Compareads

C'est ain­si que j'ai déve­lop­pé le logi­ciel Com­pa­reads durant ma thèse. Ce logi­ciel per­met de com­pa­rer deux échan­tillons méta­gé­no­mique pour iden­ti­fier et comp­ter les séquences simi­laires entre eux.
On appelle "séquences simi­laires", deux séquences (une venant de chaque échan­tillon) qui ont en com­mun au moins t k-mers. Dans la pra­tique, on uti­lise sou­vent t=2 et k=33 pour des séquences de 100 paires de bases. Donc, si une séquence de A et une séquence de B ont au moins deux sous séquences dis­tinctes de 33 nucléo­tides en com­mun, on consi­dère que ces deux séquences sont simi­laires (figure 1).

Un algo­rithme sim­pliste pour expli­quer Com­pa­reads est :

Pour  chaque  sequence  s  dans  A

{

    Pour  chaque  k-mer  m  de  s

    {

        Si  m  est  dans  B

        {

            nbkm++

        }

    }

    Si  nbkm  >=  t { afficher  s }

}



Cet algorithme permet d'identifier les séquences de A similaires à au moins une séquence de B. On appelle cette opération une demi-intersection, notée . Ce symbole signifie à la fois que c'est une demi-intersection (donc dirigée dans le sens de la flèche) et la flèche en zigzag est là pour exprimer le fait que ce résultat n'est pas exact mais légèrement sur-évalué (une sombre histoire de faux-positifs sur laquelle je reviendrai plus tard).
La figure 2 reprend cet algorithme visuellement : chaque barre horizontale noire est une séquence et les boîtes de couleurs sont les k-mers identiques entre A et B. Trouver toutes les séquences de l'ensemble A similaires à au moins une séquence de l'ensemble B peut se représenter ainsi :

Il faut noter que cette demi-intersection ne renseigne pas sur les séquences de B similaires à des séquences de A, seulement sur les séquences de A similaires à des séquences de B.
Le processus complet pour identifier toutes les séquences similaires entre A et B, donc trouver (AB) mais aussi (BA) est composé de trois demi-intersections et non deux comme on pourrait s'y attendre (figure 3). Le premier résultat (AB)* n'est pas tout à fait identique à l'image 2, encore la faute à des faux-positifs. Ce résultat est temporaire c'est pourquoi il est suivi d'un astérisque.

La première demi-intersection n'est pas exactement comme montré figure 2. Souvenez-vous qu'on indexe tous les k-mers de B (et seulement les k-mers de B) dans le filtre de Bloom. Ce faisant, pour B, on perd l'information "ce k-mers provient de cette séquence". Dès lors, on peut avoir un soucis...

Si on compare le résultat de la demi-intersection entre A et B sur la figure 2 et 4, on voit qu'il y a deux séquences de plus identifiées sur le figure 4. En effet, en utilisant exclusivement l’information "ce k-mer est présent ou absent dans B", on peut se retrouver avec une séquence de A partageant des k-mers avec plusieurs séquences distinctes de B. Autrement dit, cette séquence de A a bien t k-mers de partagés avec B, seulement ces k-mers ne sont pas sur la même séquence de B. La séquence de A est, à tord, identifiée comme similaire à une séquence de B.
Par contre, toutes les séquences de A effectivement similaires à des séquences de B sont bien présentes dans le résultat. On peut donc utiliser ce sous-ensemble de A lors de la seconde étape de Compareads (trouver les séquences de B similaires à des séquences de A) en lieu et place de A au complet. Moins de séquences à indexer signifie un gain de temps et moins de risque d'avoir des erreurs comme celle que je viens de décrire.

Finalement, on recommence le première étape, trouver les séquences de A similaires à des séquences de B, mais en utilisant le résultat de la seconde étape.

Ceci nous permet d'avoir deux résultats (AB et BA) qui utilisent eux-même le résultat d'une comparaison. Ce mécanisme un peu compliqué à comprendre à première vue permet de diminuer le nombre de séquences faussement identifiées comme similaires (faux positifs) : vous pouvez voir figure 6 qu'il y a une séquence mal identifiée en moins comparé à la figure 4. Malgré cela, il peut rester des faux positifs, comme montré  figure 6. Le taux de faux positif de Compareads est tout de même très faible et estimé à moins de 1%.
Trois scores sont alors calculés. Pour chaque échantillon, on a le nombre de séquences de cet échantillon similaires à des séquences de l’autre échantillon. Le troisième score ressemble à l'indice de Jaccard et est calculé en additionnant le nombre de séquences de A similaire à des séquences de B et le nombre de séquences de B similaire à des séquences de A, puis en divisant ce nombre par le nombre total de séquence de A et B.
Commet
Commet est l'évolution de Compareads. Il me faudrait un article entier pour l'expliquer en détail, je vais donc seulement le survoler.
Il reprend l'algorithme de Compareads mais ne calcule pas simplement la similarité entre deux échantillons : il calcule toutes les intersections entre N échantillons pris deux à deux, en optimisant l'utilisation de mémoire vive, le temps de calcul (au moins 2 fois plus rapide) et divise par 100 l'impact des sorties sur disque dur.

Commet génère plusieurs sorties graphiques, comme la "heatmap" figure 7. Sur cette carte, le vert représente une similarité entre deux échantillon faible et le rouge une similarité plus forte. On voit alors facilement quels sont les échantillons qui sont plus proche d'un point de vue séquences.
De plus, Commet autorise l'utilisation d'opérations binaires entre les résultats d'intersection de deux échantillons. Sur la figure 8, on voit en bleu les séquences de A ET NON(AB), en vert (AB) ET (AC), en orange, (BA) OU (BC) et en rouge, (CA) ET NON (C B).

Cela permet de lancer une seule fois Commet sur tous les échantillons d'un projet métagénomique (par exemple, les échantillons de l'expédition Tara Oceans), pour avoir les pourcentages de similarité entre tous les échantillons pris deux à deux (et les séquences similaires qui y correspondent). Alors, on peut presque instantanément (quelques secondes) répondre à des questions telles : quel est le pourcentage de similarité entre les séquences de cet échantillon de l'océan Arctique et les séquences de ces trois échantillons de l'océan Antarctique qui ne sont pas similaire aux séquences de ce deux échantillons des Îles Malouines ? Ou encore, quelles sont les séquences communes à la mer Rouge, l'océan Indien et l'océan Pacifique sud mais qu'on ne retrouve pas dans la Méditerranée ?
Références :
N Maillet, C Lemaitre, R Chikhi, D Lavenier, P Peterlongo. Compareads: comparing huge metagenomic experiments. BMC Bioinformatics 2012 13(Suppl 19):S10.
N Maillet, G Collet, T Vannier, D Lavenier, P Peterlongo. COMMET: comparing and combining multiple metagenomic datasets, IEEE International Conference on Bioinformatics and Biomedicine (BIBM) 2014. Télécharger l'outil.
 
Si le sujet vous intéresse, sachez qu'un nouvel algorithme nommé SimKa est actuellement en construction à INRIA Rennes par Gaëtan Benoît pour estimer la similarité entre N jeux de données métagénomique encore plus rapidement. Il a gagné le prix du meilleur poster de JOBIM2015.
Un grand merci à Bérénice Batut et m4rsu pour avoir relu et corrigé cet article :)
      
          
        Vous avez aimé ? Dites-le nous !      
            
      
                  
                      
                  
                      
                  
                      
                  
                      
                  
                      
              
    
    
          Soumettre        
      Moyenne : 0 / 5. Nb de votes : 0    
    
      Pas encore de vote pour cet article.