Métagénomique et Graphique de Recrutement : chercher l'aiguille dans la botte de foin !

Le massacre des Albigeois
La croi­sade contre les Albi­geois, vers 1208. Manus­crit tiré des "Chro­nique de France ou de Saint Denis", 1325–1350. Source : Wiki­me­dia Com­mons.

Un inté­rêt cer­tain se déve­loppe pour la com­pré­hen­sion d'écosystèmes com­plexes comme les eaux, les sols ou encore les micro­biomes (plus de 90% des cel­lules du corps humain sont en fait les bac­té­ries qui peuplent son tube diges­tif). Plu­sieurs dif­fi­cul­tés se posent :

  1. la plu­part des orga­nismes pré­sents ne sont pas faci­le­ment culti­vables et donc indi­vi­dua­li­sables.
  2. tous les orga­nismes ne sont pas pré­sents en quan­ti­té égale et, par exten­sion, d'importance égale dans le fonc­tion­ne­ment de l'écosystème, sans cor­ré­la­tion stricte entre abon­dance et impor­tance.
  3. cer­tains de ces éco­sys­tèmes sont extrê­me­ment riches et des cen­taines d'espèces peuvent se cacher dans quelques grammes de ter­reau ou quelques mil­li­litres d'océan.

Face à une tâche aus­si ingrate que pas­sion­nante, les cher­cheurs ont trou­vé une solu­tion dans la célèbre phrase asso­ciée à la Croi­sade des Albi­geois : « Tuez-les tous, Dieu recon­naî­tra les siens », en l'adaptant légè­re­ment :

« Séquen­cez-les tous, un doc­to­rant fera le tri »

Ain­si naquît la Méta­gé­no­mique et de grands pro­jets comme le séquen­çage du rumen des vaches ou du micro­biome humain, le pro­jet METASOIL, le Glo­bal Ocean Sur­vey et plus récem­ment le pro­jet TARA.

Il convient de dis­tin­guer deux pos­si­bi­li­tés en méta­gé­no­mique :

  1. séquen­cer un mar­queur géné­tique capable d'identifier les orga­nismes pré­sents à l'échelle du genre voire de l'espèce comme le mar­queur V9 (cor­res­pon­dant à l'ADNr 18S pour les euca­ryotes, 16S pour les pro­ca­ryotes)
  2. abu­ser des tech­niques modernes de séquen­çage très haut débit : extraire la tota­li­té des acides nucléiques de l'échantillon, tout décou­per, tout séquen­cer. Pour l'analogie, ima­gi­nons une virée à la Biblio­thèque Natio­nale, un décou­page inten­sif de l'ensemble des pages de cha­cun des livres (dont cer­tains sont ‑au grand bon­heur des nom­breux fans d'Har­ry Pot­ter- en plu­sieurs copies), la réa­li­sa­tion d'un gros tas rela­ti­ve­ment stable, et pour finir se lan­cer, avec espoir, témé­ri­té et un peu de ruban adhé­sif, dans la recons­ti­tu­tion de tous ces chefs d’œuvre (si si, même ceux qu'on n'aime pas, comme 50 Shades of Grey ou la saga Twillight).

Fort heu­reu­se­ment, rares sont les cher­cheurs qui osent se lan­cer en quête d'une com­pré­hen­sion exhaus­tive des popu­la­tions d'un éco­sys­tème com­plexe. Si la quête se révèle rela­ti­ve­ment réa­li­sable lorsqu'on ne séquence qu'un mar­queur uni­ver­sel et qu'il s'agit d'identifier "qui est là", ten­ter de recons­ti­tuer les génomes via une étape d'assemblage se révèle un défi. Bien sou­vent et c'est ce que nous allons abor­der ici, le cher­cheur lamb­da ne s'intéressera qu'à détec­ter l'abondance de son picoeu­ca­ryote ado­ré ou de ses bac­té­ries favo­rites. Pour cela un outil astu­cieux et très pra­tique existe : le Gra­phique de Recru­te­ment (libre­ment tra­duit du Recruit­ment Plot que vous trou­ve­rez dans la lit­té­ra­ture)

L'idée étant sim­ple­ment de recru­ter ‑par ali­gne­ment- des reads de séquen­çage du méta­gé­nome sur un ou plu­sieurs génomes de réfé­rence. Et de dis­tin­guer ces reads par leur pour­cen­tage d'identité avec la séquence cible. Lorsqu'une quan­ti­té signi­fi­ca­tive de reads s'alignent sur le génome d'une souche de réfé­rence, et ce avec un pour­cen­tage d'identité suf­fi­sam­ment impor­tant, on peut émettre l'hypothèse que cet orga­nisme ou une varié­té d'organismes très proches sont pré­sent dans l'échantillon consi­dé­ré.

Ingré­dients de la recette du jour : Le Gra­phique de Recru­te­ment

  1. Quelques dizaines de mil­lions de "reads" d'une taille hon­nête (pas moins de 100 paires de bases), au mieux une banque Illu­mi­na "Pair-End" de 2 x 100 pb che­vau­chants voire de longs reads 454 ou Ion­Tor­rent.
  2.  Quelques génomes de réfé­rence, repré­sen­ta­tifs de la diver­si­té du groupe dont vous sou­hai­tez détec­ter la pré­sence
  3.  Un gros, gros ordi­na­teur (ou bien mieux un clus­ter de cal­cul)
  4. Quelques com­pé­tences en "lan­çage de BLAST" et dans un lan­gage de pro­gram­ma­tion "User Friend­ly" (type Bash, Python ou Perl) ain­si qu'un peu de R pour faire de beaux gra­phiques, mais Excel reste uti­li­sable pour les fri­leux
  5. Beau­coup d'enthousiasme !

Dif­fi­cul­té : 3 /​ 5     Rien d'extravagant ici, mais néces­site une cer­taine sou­plesse devant une ligne de com­mande

Méthode en 3 étapes :

1/​ Ali­gner l'ensemble des reads du méta­gé­nome sur l'ensemble des génomes "de réfé­rence" à l'aide de BLAST, en auto­ri­sant des "hits" aus­si faibles que 55% d'identité

L'idée étant de trou­ver les para­mètres appro­priés :

blastall -p blastn -i myReads.fasta -d myDatabase -o outfile.tab -m 8 -G 8 -E 6 -r 5 -q -4 -W 8 -e 1 -F "m L" -U T

Avec ‑W 8 : la taille du "mot", ‑G 8 ‑E 6 ‑r 5 ‑q 4 : les coûts d'ouverture et d'élongation de gap, la récom­pense pour un match et la péna­li­té pour un mis­match, -e 1 -F "m L" -U T : quelques filtres pour allé­ger les fichiers de sor­tie, et non moins impor­tant ‑m 8 pour une sor­tie tabu­lée facile à mani­pu­ler. Ces para­mètres ne sont qu'une sug­ges­tion, il convien­dra de les adap­ter à ses propres don­nées.

2/​ Net­toyer les résul­tats (à l'aide de Bash, Python ou Perl, ou autre…)

Il s'agit de conser­ver les hits BLAST qui vous semblent per­ti­nents, non pas en termes de pour­cen­tage d'identité, mais en terme de lon­gueur de l'alignement. J'essaye par exemple d'avoir au moins 80% de la séquence de mon read qui soit effec­ti­ve­ment ali­gnée. À ce stade, vous vou­lez savoir si un orga­nisme de réfé­rence don­né est pré­sent dans le méta­gé­nome. Il ne s'agit pas de faire de com­pé­ti­tion entre vos orga­nismes de réfé­rence, je vous sug­gère donc de conser­ver éga­le­ment l'information de reads qui s'alignent sur plu­sieurs de vos génomes (en l'occurence, consi­dé­rant que vos génomes sont proches, cela devrait se pro­duire régu­liè­re­ment). Par contre, n'hésitez pas à ne conser­ver que le meilleur hit de chaque read pour cha­cune de  vos souches de réfé­rence.

Il s'agit ensuite d'extraire trois infor­ma­tions clefs : le génome de réfé­rence où le read s'est ali­gné, le pour­cen­tage d'identité du hit et la loca­li­sa­tion de ce hit sur le génome de réfé­rence (pas de pitié, faites la moyenne du Stop et du Start).

3/​ Construc­tion du Gra­phique de Recru­te­ment

Il s'agit enfin, pour un orga­nisme de réfé­rence don­né de géné­rer le gra­phique sui­vant : en abs­cisse la posi­tion sur le chro­mo­some de votre orga­nisme (s'il y en a plu­sieurs, il fau­dra faire plu­sieurs graphes ou bien les conca­té­ner), et en ordon­née, pour cha­cun des reads son pour­cen­tage d'identité. Petit conseil : il n'est pas utile d'afficher 2 mil­lions de points sur un gra­phiques, sélec­tion­nez tout au plus 100 000 points à affi­cher.

Résul­tats :

Gra­phique de Recru­te­ment pour la Souche A (G.Farrant, Tous droits réser­vés)

Gra­phique de Recru­te­ment pour la Souche B (G.Farrant, Tous droits réser­vés)

Dans ces exemples, j'ai fusion­né les pro­fils de recru­te­ment des méta­gé­nomes de trois sta­tions de pré­lè­ve­ment sur les génomes de deux souches bac­té­riennes A et B. Cette com­bi­nai­son de 3 condi­tions, asso­ciées à des cou­leurs dif­fé­rentes faci­lite la com­pa­rai­son des sta­tions. L'histogramme à droite du gra­phique repré­sente l'abondance de points sur des tranches de pour­cen­tage d'identité.

On observe que qua­si­ment aucun read ne s’aligne sur la souche A avec une iden­ti­té supé­rieure à 90% alors que sur la souche B, une quan­ti­té signi­fi­ca­tive de reads s’aligne avec plus de 90% d’identité, ce qui per­met de conclure que des orga­nismes géné­ti­que­ment proches de la souche B sont pré­sents dans l’échantillon étu­dié. Les points recru­tés entre 55 et 80% d'identité cor­res­pon­draient à du recru­te­ment non-spé­ci­fique.

On observe éga­le­ment quelques arte­facts poten­tiel­le­ment inté­res­sants tels des sites où le signal est très riche, l'intuition amène à pen­ser qu'il s'agit de régions conte­nant des gènes de ménage lar­ge­ment conser­vés entre orga­nismes (type ARN ribo­so­miques) ain­si que des régions "vides", dont une très mar­quée pour la Souche 2, qui doivent être inter­pré­tés (îlots géno­miques ?).

Le Gra­phique de Recru­te­ment est un pre­mier outil d'analyse de la pré­sence d'organismes de réfé­rence dans un méta­gé­nome. Il per­met entre autres de défi­nir un seuil de pour­cen­tage d'identité prompt à dis­cri­mi­ner les familles d'organismes. Il per­met éga­le­ment de révé­ler des par­ti­cu­la­ri­tés riches d'interprétations. Pour plus d'infos et du sub­ti­li­tés, je vous sug­gère la publi­ca­tion du Glo­ba­lO­cean­Sur­vey (GOS) qui s'octroie la pater­ni­té de la méthode.

Je tiens à remer­cier Mali­cia, Haut­bit, Yoann M., Nico M., Mica et Guillaume Col­let pour leurs com­men­taires avi­sés.




Commentaires

7 réponses à “Métagénomique et Graphique de Recrutement : chercher l'aiguille dans la botte de foin !”

  1. Pour mini­mi­ser le net­toyage, j'aurai effec­tuer un lan­çage de blast avec une valeur de ‑e bien infé­rieure à 1 😛

    Sinon, les sor­ties gra­phiques des graphes de recru­te­ment sont de toute beau­té. Pos­si­bi­li­té de par­tage du code ?

    1. Oui, je suis assez d'accord sur le prin­cipe, mais en véri­té, je ne vou­lais pas prendre le risque de lou­per quelque chose et n'étant pas limi­té par l'espace disque, j'ai pré­fé­ré faire du net­toyage post-BLAST.

      En ce qui concerne les sor­ties gra­phiques, je ne suis pas en mesure de four­nir le code, mais il n'y a clai­re­ment rien d'ésotérique vu que j'ai pu le faire tout seul.

  2. Cette repré­sen­ta­tion est vrai­ment très inté­res­sante, elle pour­rait per­mettre d'obtenir un visuel com­pa­ra­tif plus com­plet des génomes connus dans ce type de jeux de don­nées qui sont quand même assez consé­quents et "fouillis" à appré­hen­der de prime abord.

    Quelques ques­tions (je n'ai pas encore lu la publi donc par­donne-moi si ça y est men­tion­né) :
    — Lorsqu'on observe les 3 pro­fils de 3 sites dif­fé­rents sur un même génome et un même graphe, inclues-tu une étape de nor­ma­li­sa­tion de la quan­ti­té des reads ?
    — 90% d'identité pour vali­der une souche me parait peu, com­ment as-tu fixé ton seuil ? (je serais plus encline à ne vali­der que l'espèce avec un tel pour­cen­tage)

    Je sens que je pour­rais m'inspirer de tout ceci pour une repré­sen­ta­tion com­plé­men­taire à cer­taines de mes ana­lyses, mer­ci en tout cas de l'avoir mise en avant, il va fal­loir que je creuse !

    PS : Petite note de pinailleuse, il serait judi­cieux de mettre à jour la ligne de com­mande blast avec la nou­velle nomen­cla­ture des der­nières ver­sions (blastn ‑que­ry, ‑db, ‑out, ‑word_​size, ‑eva­lue, …) (je sais, pas­ser de l'un à l'autre est une gym­nas­tique bien pénible, je me plante encore sou­vent…)

    1. Tout d'abord, mer­ci pour ton inté­rêt.

      Alors non, il n'y pas d'étape de nor­ma­li­sa­tion pour la concep­tion des Gra­phiques de Recru­te­ment car je me contente d'afficher les 100 000 pre­miers 'hits' de chaque sta­tions. J'ai pris le temps de voir les résul­tat en affi­chant la tota­li­té de mes points, c'est pas­sa­ble­ment illi­sible, et en échan­tillon­nant les 100 000 points de façon aléa­toire, le gra­phique est rela­ti­ve­ment stable à l’échantillonnage.
      La com­pa­rai­son fine des sta­tions est dif­fi­cile sur le Gra­phique de recru­te­ment, per­son­nel­le­ment, je m'en sers plus comme d'une base pour détec­ter les sus-nom­més arte­facts et com­pa­rer leur pré­sence dans les dif­fé­rentes sta­tions.
      Com­pa­rer les sta­tions de façon quan­ti­ta­tive est dif­fi­cile car on ne connaît pas sys­té­ma­ti­que­ment le volume ou la masse de l'échantillon de méta­gé­nome. Par contre, on peut ‑et on doit- faire du semi-quan­ti­ta­tif et com­pa­rer les pro­fils de dis­tri­bu­tion entre les sta­tions. Typi­que­ment, on peut faire du recru­te­ment sur les génomes de 36 souches appar­te­nant à 17 clades dif­fé­rents, il est pos­sible, en dénom­brant les meilleurs hits BLAST pour cha­cun des clades d'en esti­mer la pro­por­tion dans le méta­gé­nome, et ce pro­fil semi-quan­ti­ta­tif est quand à lui com­pa­rable entre les sta­tions.

      Concer­nant la ques­tion du seuil, il n'y pas de bonne réponse. Une stra­té­gie pos­sible consiste à éta­blir arbi­trai­re­ment un seuil, puis extraire les séquences cor­res­pon­dantes, les 'blas­ter' sur nr et voir le ratio de séquences qui ont un 'best hit' sur l'organisme de réfé­rence, adap­ter le seuil, et recom­men­cer, mais nous sommes d'accord ‑je l'espère- pour pen­ser que LA souche de réfé­rence n'existe que dans la lit­té­ra­ture et que le milieu étu­dié contient au mieux un spectre rela­ti­ve­ment hété­ro­gène d'organismes proches. Fina­le­ment on en revient à la ques­tion de défi­ni­tion d'espèce.

      Ravi que ça t'ai ins­pi­ré. (Pour le pinaillage, je dirais "You're rai­sing an inter­es­ting concern" ^^)

  3. Avatar de AnoukLy
    AnoukLy

    Bon­jour,

    Tout d'abord mer­ci beau­coup pour cette article ain­si que pour l'ensemble de votre site qui est très inté­res­sant.

    Je me per­mets de sol­li­ci­ter votre aide car j'aimerais effec­tuer une mani­pu­la­tion simi­laire sur mes échan­tillons, mais hélas je débute en méta­gé­no­mique.
    J'aurais vou­lu savoir par quel moyen vous déci­diez des valeurs pour les dif­fé­rents para­mètres du blast. Ensuite auriez vous plus de ren­sei­gne­ments sur la manière d'effectuer le net­toyage sur les résul­tats du blast ? J'ai essayé de trou­ver de la docu­men­ta­tion ou des tuto­riels, mais je ne trouve rien d'accessible ou de très pré­cis.

    Par avance mer­ci beau­coup pour votre aide !

    1. Avatar de NiGoPol
      NiGoPol

      Bon­jour,
      Déso­lé de n'avoir pris le temps de répondre plus tôt.
      Ravi que cela vous plaise =)

      Bien, alors concer­nant BLAST, je pense qu'il fau­drait édi­ter un erra­tum. En fait, à l'époque du début de ma thèse, j'utilisais BLASTALL (une ver­sion plus ancienne, mais pas obso­lète de BLAST. Avec les para­mètres par défaut, BLASTALL et la der­nière ver­sion du logi­ciel (appe­lée "BLAST+") se com­porte comme MEGABLAST et génère des ali­gne­ments géné­ti­que­ment proches com­pris glo­ba­le­ment entre 80 et 100% d'identité, et donc inuti­li­sables pour la construc­tion des­dits "recruit­ment plots", ces der­niers se basent sur des ali­gne­ments entre 50 et 100%.
      Les para­mètres "sen­sibles" de blast pré­sen­tés dans cet article per­mettent de for­cer MEGABLAST à géné­rer des ali­gne­ments dis­tants jusqu'à ~50% d'identité.

      En fait, j'ignorais à l'époque de la rédac­tion l'existence d'une option de BLAST+: ‑task blastn qui lui per­met de se com­por­ter comme un BLASTN nor­mal. L'utilisation des options "sen­sibles" de BLAST per­mettent de géné­rer plus d'alignements.

      J'essaierai de faire un article plus détaillé sur la ques­tion.

      1. Avatar de AnoukLy
        AnoukLy

        Mer­ci pour ces infor­ma­tions com­plé­men­taires !

Laisser un commentaire