Soirée BED & FASTA !

Après la petite his­toire de l’analyse des séquences d’ADN, voi­ci un tuto­riel pour apprendre quelques trucs et astuces dans ce domaine.

Bio­lo­giste en mal de connais­sances de pro­gram­ma­tion ou pro de R, vous trou­ve­rez ici de quoi vous amu­ser avec un fichier Fas­ta ou un Bed.

Nous allons voir com­ment faire un ali­gne­ment mul­tiple de séquences (voir ici pour un état de l'art des logi­ciels exis­tants), pro­duire un logo, ana­ly­ser les motifs enri­chis dans des pics de CHiP, et plus encore. Vous aurez la pos­si­bi­li­té de tra­vailler soit avec des inter­faces web, soit avec des lignes de codes R et dif­fé­rentes librai­ries Bio­con­duc­tor (que nous vous avions déjà pré­sen­té).

Vous pou­vez télé­char­ger l'archive conte­nant les fichiers néces­saire pour le tuto­riel ici.

Vous avez dit Fasta ?

Un Fas­ta est un fichier texte conte­nant une ou plu­sieurs séquences d’ADN, ARN ou acides ami­nés. Chaque séquence se pré­sente sous la forme sui­vante :

Le signe ‘>’ repré­sente la ligne avec le nom de la séquence. Il contient en géné­ral un iden­ti­fiant ou la loca­li­sa­tion géno­mique de la séquence. Ce for­mat est très utile car il est com­pa­tible avec beau­coup d’outils dis­po­nibles sur inter­net (inter­faces web) et de nom­breuses librai­ries de fonc­tions, notam­ment adap­tées à R, Python et Perl. Vous pou­vez télé­char­ger des séquences en for­mat Fas­ta depuis Gene­Bank et UCSC, entre autres.

Le for­mat Bed est aus­si un fichier texte qui contient plu­sieurs loca­li­sa­tions sur le génome. A la base, un Bed contient trois colonnes : le nom du chro­mo­some, la posi­tion de départ de la région et la posi­tion de fin de la région.  Dans cer­tains cas le fichier Bed contient un nom pour chaque région, un score (i.e. une quan­ti­fi­ca­tion) ain­si que le brin posi­tif ou néga­tif pour l’orientation géno­mique de la région.

Vous pou­vez trou­ver d’autres infor­ma­tions sur les for­mats de don­nées géno­miques sur le site UCSC.

Donc main­te­nant que vous avez les bases, pas­sons aux choses sérieuses !

Identification de la séquence suspecte

Allez sur le site de BLAST, et cli­quez sur "nucleo­tide blast". Copiez et col­lez La séquence Fas­ta ci-des­sus pour iden­ti­fier la séquence sans la ligne avec le signe ’>’.

Alors ? Vous avez trou­vé ?

Effec­ti­ve­ment il s’agit d’une par­tie de la séquence du récep­teur aux glu­co­cor­ti­coïdes.

Cette hor­mone est un outil thé­ra­peu­tique impor­tant dans le cadre des mala­dies auto-inflam­ma­toires et dans l’inhibition de la crois­sance cel­lu­laire de nom­breux can­cers. Les glu­co­cor­ti­coïdes sont syn­thé­ti­sés et sécré­tés par le sys­tème ner­veux cen­tral  et les glandes cor­ti­co-sur­ré­nales, sui­vant un rythme cir­ca­dien. Ils sont impli­qués dans la régu­la­tion de nom­breux pro­ces­sus bio­lo­giques, entre autres la glu­co­néo­ge­nèse, la gly­co­lyse, le méta­bo­lisme des acides gras et la réponse immu­ni­taire et inflam­ma­toire.

Les glu­co­cor­ti­coïdes se fixent aux récep­teurs des glu­co­cor­ti­coïdes (GR) dans le cyto­plasme de la cel­lule. Le com­plexe récep­teur-ligand for­mé pénètre ensuite dans le noyau cel­lu­laire où il se fixe à de nom­breux élé­ments de réponse aux glu­co­cor­ti­coïdes dans la région du pro­mo­teur des gènes-cibles. Le récep­teur, ain­si fixé à la molé­cule d'ADN, inter­agit avec les fac­teurs de trans­crip­tion basiques, pro­vo­quant une aug­men­ta­tion de l'expression de gènes-cibles spé­ci­fiques. Ce pro­ces­sus est appe­lé tran­sac­ti­va­tion. Le méca­nisme inverse est appe­lé trans­ré­pres­sion. Le récep­teur hor­mo­nal acti­vé inter­agit avec des fac­teurs de trans­crip­tion spé­ci­fiques et pré­vient la trans­crip­tion des gènes-cibles. Les glu­co­cor­ti­coïdes ont une action pléio­tro­pique : ils activent cer­tains gènes (tran­sac­ti­va­tion) et en répriment d’autres (trans­ré­pres­sion) en fonc­tion de la struc­ture des pro­mo­teurs ciblés et de leur action coor­don­née avec des cofac­teurs.

Dans un autre contexte, cer­tains bio­lo­gistes uti­lisent blastn pour savoir si un insert s’est cor­rec­te­ment  intro­duit dans un plas­mide lors d’un clo­nage. L’on effec­tue cette véri­fi­ca­tion par élec­tro­pho­rèse sur gel d’agarose, puis ampli­fi­ca­tion PCR et enfin par séquen­çage. Le Fas­ta issu du séquen­çage peut ensuite être ana­ly­sé avec BLAST pour retrou­ver la posi­tion exacte de l’insert sur le plas­mide pro­duit.

Alignement multiple de sites GR et Logo

Ren­dez-vous main­te­nant sur le site de Clus­tal Ome­ga pour ten­ter un ali­gne­ment mul­tiple de 17 séquences issues de mRNA de GR de dif­fé­rentes espèces et connaître la dis­tance phy­lo­gé­né­tique entre elles. Pour cela, il vous suf­fit de char­ger le fichier Fas­ta "GR_phylo.fasta" , de pré­ci­ser qu’il s’agit d’ADN et de spé­ci­fier le for­mat de sor­tie en Fas­ta. Vous ver­rez rapi­de­ment le résul­tat de l’alignement (ci-des­sous, en for­mat Clus­tal pour un petit inter­valle).

En cli­quant sur « phy­lo­ge­ne­tic tree » vous devriez avoir l’image sui­vante :

L’arbre per­met des savoir quelles sont les séquences les plus proches. En uti­li­sant ce type d’analyse, on peut com­pa­rer dif­fé­rentes espèces qui par­tagent cer­tains gènes et ain­si retra­cer leur « his­toire évo­lu­tive ». Par ailleurs, si vous dési­rez vous attar­der un peu plus sur la construc­tion de votre arbre et sou­hai­tez une plus forte robus­tesse, je vous invite alors à consul­ter notre article sur la construc­tion des arbres phy­lo­gé­né­tiques.

On peut aus­si uti­li­ser l’alignement pour créer un logo avec webLo­go pour savoir quelles sont les bases « conser­vées ». Voi­ci le résul­tat pour la région 3500–3700 de l’alignement :

Dans cette ana­lyse, nous pou­vons obser­ver 3 groupes de séquences GR d’espèces proches, au sens de la  dis­tance géné­tique. Le logo nous per­met de consta­ter la conser­va­tion de cer­taines régions de l’alignement mul­tiple. Le but de cet exemple est de com­prendre com­ment fonc­tionnent ces dif­fé­rentes inter­faces web.

Analyse de motifs enrichis dans les pics de CHiP-seq du récepteur aux glucocorticoïdes (GR)

En géné­ral, lors d’analyses sur des don­nées CHiP, nous avons plus de 1000 sites. L’analyse des séquences est alors plus com­pli­quée car les ali­gne­ments mul­tiples sont trop coû­teux en matière de cal­culs. D’autres algo­rithmes ont vu le jour pour répondre plus spé­ci­fi­que­ment à cette ques­tion. Il s’agit d’algorithmes de recherche de motifs de novo  (i.e. MEME) ou encore de scans de séquences à par­tir de PWM (i.e. FIMO).

Nous allons pro­cé­der à une ana­lyse de don­nés de 1000 pics de CHiP-seq de GR pro­ve­nant d’un papier de G.Hager. Télé­char­gez le fichier  "GSE46047_GR_1000peaks_mouse_liver_mm9.txt" conte­nant la loca­li­sa­tion des pics (for­mat .bed). (voir archive)

Allez sur UCSC et sélec­tion­nez le génome de sou­ris mm9. Cli­quez ensuite sur le bou­ton « add cus­tom track » et char­gez le fichier de loca­li­sa­tion des pics. Allez ensuite sur l’onglet « Tables » pour conver­tir le fichier en Fas­ta.

Sau­vez ensuite le fichier Fas­ta sous "GR1000.fa", puis allez sur MEME-CHiP pour sou­mettre le Fas­ta à la recherche de motifs de novo. Uti­li­sez les para­mètres sui­vants :

L’exécution du pro­gramme MEME peut être assez longue, mais vous rece­vrez par e‑mail le lien qui vous fera accé­der à votre résul­tat. MEME-CHIP vous four­ni­ra les motifs enri­chis et leur dis­tri­bu­tion par rap­port au centre des pics. La docu­men­ta­tion de MEME-CHIP contient beau­coup d’informations per­met­tant de com­prendre les résul­tats de l’analyse. Nous allons voir main­te­nant com­ment faire à peu près la même ana­lyse en uti­li­sant R.

Un peu de code R pour l’analyse de séquences ADN

Dans un pre­mier temps il faut ins­tal­ler dif­fé­rents packages R :

Le package BSgenome.Mmusculus.UCSC.mm9 fait 750MB, donc l’installation risque de prendre un peu de temps ! Le café s’impose à ce stade…

Puis char­gez les biblio­thèques :

Il faut ensuite char­ger les don­nées des pics de CHiP-seq de GR :

Les com­mandes sui­vantes per­mettent de créer 2 fichiers Fasta(20 sites et 1000 sites) à par­tir des pics (fichier BED, si jamais vous n'avez pas de fichier FASTA déjà à dis­po­si­tion).

Nous allons faire un ali­gne­ment mul­tiple sur le fichier Fas­ta conte­nant 17 séquences de mRNA d’espèces dif­fé­rentes de GR grâce à “muscle”.

En uti­li­sant l’alignement,  nous pou­vons faire une ana­lyse de la dis­tance en terme de simi­la­ri­té entre les séquences.

Le résul­tat de cette ana­lyse peut être repré­sen­té sous forme de heat­map.

Un den­dro­gramme peut être obte­nu comme suit :

Le den­dro­gramme ain­si que la heat­map nous montrent que les séquences des 17 mRNA de GR ont une simi­la­ri­té rela­ti­ve­ment impor­tante. On peut dis­cer­ner 3 groupes, avec res­pec­ti­ve­ment des séquences pro­ve­nant d’espèces de pois­sons, de mam­mi­fère et un troi­sième groupe qui com­porte des séquences mRNA de GR pro­ve­nant de rat, de poule et de gre­nouille.

Le logo de l’alignement est obte­nu comme suit.

Le but de cet exemple est sim­ple­ment de com­prendre com­ment repré­sen­ter un logo avec le package motif­Stack à par­tir d’un objet DNAS­tring­Set pro­ve­nant du package Bios­trings.

Une analyse plus réaliste des données CHiP-seq de GR

L’analyse de séquences  de CHiP GR implique, la recherche de motifs de novo. Cette ana­lyse peut se faire sur un nombre plus impor­tant de sites. De plus elle per­met de trou­ver des motifs de cofac­teurs poten­tiels  proches du site de liai­son de la pro­téine GR. Cette ana­lyse est pos­sible grâce au packages rGA­DEM et MotIV. L’analyse est assez longue, donc vous pou­vez lan­cer les com­mandes sui­vantes et sor­tir boire une bière avec des amis !

Vous devriez obte­nir l’image ci-des­sous. Le motif le plus enri­chi est NR3C1, qui cor­res­pond aux sites de liai­sons de la pro­téine GR. Nous trou­vons d’autres motifs comme ESR2, motif de liai­son des récep­teurs à œstro­gènes, et STAT3 et PPARG:RXR qui cor­res­pondent aux sites de liai­son de cofac­teurs connus de GR.

L’analyse de sites GR de CHiP montre un lien avec Stat3. De ce fait, nous pou­vons ana­ly­ser la dis­tri­bu­tion de sites de Stat3 et de GR autour du centre des pics de CHiP-seq grâce aux com­mandes sui­vantes :

Enfin nous pou­vons voir com­bien de sites de liai­son de GR che­vauchent des sites de liai­son de Stat3, et ana­ly­ser ain­si la dis­tance rela­tive entre GR et Stat3.

C’est tout pour aujourd’hui ! Grâce à ces quelques lignes de codes et ces outils dis­po­nibles sur le web, vous pour­rez ana­ly­ser de manière spé­ci­fique et sys­té­ma­tique vos don­nées CHiP-seq et vos séquences.

Encore une petite remarque :

L’algorithme de BLAST est en géné­ral appe­lé en ligne de com­mande depuis R, à condi­tion d’avoir l’outil ins­tal­lé et acces­sible depuis le ter­mi­nal. Je n’ai pas trou­vé de package pou­vant l’exécuter autre­ment, donc n’hésitez pas à faire des sug­ges­tions (BoS­SA n’a pas fonc­tion­né). Vous pou­vez invo­quer BLAST depuis R avec les com­mandes ci-des­sous (non-tes­tées).

References :

1) Chro­ma­tin acces­si­bi­li­ty pre-deter­mines glu­co­cor­ti­coid recep­tor bin­ding pat­terns, Sam John, Peter J Sabo, Robert E Thur­man, Myong-Hee Sung, Simon C Bid­die, Tho­mas A John­son, Gor­don L Hager & John A Sta­ma­toyan­no­pou­los ,NATURE GENETICS, 2010

2) MEME-ChIP : motif ana­ly­sis of large DNA data­sets , Macha­nick P, Bai­ley TL.,Bioinformatics. 2011

3) GADEM : a gene­tic algo­rithm gui­ded for­ma­tion of spa­ced dyads cou­pled with an EM algo­rithm for motif dis­co­ve­ry, Li L, J Com­put Biol. 2009

4) An inte­gra­ted pipe­line for the genome-wide ana­ly­sis of trans­crip­tion fac­tor bin­ding sites from ChIP-Seq ,E Mer­cier, A Droit, L Li, G Robert­son, X Zhang, R Got­tar­do (2011). PLoS ONE. 6(2): e16432. doi:10.1371/journal.pone.0016432

5) STAMP : a web tool for explo­ring DNA-bin­ding motif simi­la­ri­ties, S. Maho­ny, P.V. Benos, Nucl Acids Res, (2007) 35:W253-258

6) DNA fami­lial bin­ding pro­files made easy : com­pa­ri­son of various motif ali­gn­ment and clus­te­ring stra­te­gies, S Maho­ny, PE Auron, PV Benos, PLoS Com­pu­ta­tio­nal Bio­lo­gy (2007) 3(3):e61

7) seq­Lo­go : Sequence logos for DNA sequence ali­gn­ments. R package Bio­con­duc­tor

8) MotIV : Motif Iden­ti­fi­ca­tion and Vali­da­tion. Eloi Mer­cier and Raphael Got­tar­do (2010). R package Bio­con­duc­tor.

9) motif­Stack : Plot sta­cked logos for single or mul­tiple DNA, RNA and ami­no acid sequence, Jian­hong Ou, Michael Brod­sky, Scot Wolfe and Lihua Julie Zhu ,R package Bio­con­duc­tor.

10) Wiki­pe­dia, article glu­co­cor­ti­coïde

Un grand mer­ci aux relec­teurs : Wocka, Yoann M, nal­lias, Julien Dela­fon­taine et Aline Jac­cot­tet.
Les images uti­li­sées pour les besoins de l'article sont toutes de l'auteur et sont dis­po­nibles en CC-by-SA 3.0.