Comment je me suis posé la question.
Chez les eucaryotes, l'ADN est organisé en domaines plus ou moins compactés, avec des taux de transcription plus ou moins élevés, et qui sont marqués différentiellement par un certain nombre de marques épigénétiques (méthylation de l'ADN, modifications post-traductionnelles des histones, variants d'histones, etc.). Il est fréquent d'essayer de corréler le niveau d'expression des gènes avec la présence ou l'absence d'une marque épigénétique à proximité des sites d’initiations de la transcription (raccourcis en TSS, pour transcription start sites). Par exemple, pour reprendre la carte de chaleur de la figure 1C de l'article présentant la suite d'outils deepTools2 :
Pour générer ce type de figures, il est nécessaire d'être en possession d'une liste de sites d’initiation de la transcription. Et donc, de se poser cette question fondamentale :
Qu'est-ce qu'un site d’initiation de la transcription ?
La définition est triviale : un site d'initiation de la transcription désigne une position génomique où le processus de transcription (synthèse d'un ARN depuis une matrice d'ADN) s'initie. En général, il est surtout fait référence aux sites d'initiations de la transcription des gènes, et/ou d'éléments transposables. Souvent l'ARN produit est coiffé au niveau de la base correspondante au site d'initiation de la transcription. Il ne faut pas confondre le site d’initiation de la transcription avec le site d'initiation de la traduction : les deux étant séparés sur les ARN messagers par la région 5' non traduite (ou 5' untranslated region, 5'UTR).
Petite subtilité : si la première base transcrite se désigne par le "+1" de transcription, la base la précédant sur le brin d'ADN est par convention appelée le "-1" de transcription, donc sans base "0" de transcription.
Comment identifier les sites d'initiations de la transcription ?
Je connais principalement deux méthodes expérimentales d'identification de ces sites. Sans rentrer dans les détails :
- la laborieuse méthode d'extension d'amorces, à faible débit.
- la plus récente méthode CAGE, pour Cap Analysis Gene Expression, un acronyme désignant un processus de capture des extrémités coiffées des ARN, suivi de leur séquençage haut débit. C'est cette méthode qui a été privilégiée par le formidable consortium Fantom5.
Les deux approches ont leurs limites. Outre la différence évidente de débit, citons le fort taux de décrochage des transcriptases inverses, facteur limitant de la méthode d'extension d'amorce, mais aussi la complexité du protocole expérimental du CAGE, ainsi que son incapacité à détecter les sites d'initiations non coiffés au niveau des ARNs.
Oui, ça OK, mais on la télécharge où la liste des TSS ?
Je reconnais bien là votre pragmatisme de bioinformaticiens. Pour nous, nul besoin de réflexions bio-philosophiques poussées : un gène, c'est ce qui est présent dans notre fichier d'annotation des gènes, et cette approche nous convient la plupart du temps ! Et les sites d'initiations de la transcription se trouvent au début des gènes.
(Ouvrons une parenthèse tout de même, pour rappeler qu'il est indispensable de bien connaître les limites actuelles de ces annotations. Ainsi, si l'humain a environ 20 000 ARN longs non-codants (lncRNA) et l’opossum seulement 8 000, ce n'est sans doute pas que nos 12 000 lncRNA de plus nous rendent supérieurement intelligents, mais plutôt que de nombreux lncRNA ne sont pas encore identifiés chez l'opossum !)
Pour les espèces sur lesquelles j'ai le plus travaillé (humain, souris), la source de référence des annotations des gènes est l'organisation GENCODE, dont les annotations sont utilisées par Ensembl.
Mais, et c'est là où je voulais en venir en commençant cet article : avec l'amélioration des annotations, les modèles de gènes comportent de plus en plus de transcrits différents. Prenons par exemple le cas du gène BRCA1 (pour BReast CAncer 1, un gène connu pour son implication dans les cancers du sein et de l'ovaire).
Où est le TSS de BRCA1 ?
Regardons à quoi ressemble les modèles de transcrits du gène BRCA1 dans la version actuelle des annotations Gencode (la 27). C'est un gène situé sur le brin moins, son début est donc sur la droite, et il est transcrit vers la gauche. D'après l'image Ensembl :
J'ai dénombré à l’œil 30 transcrits différents, auxquels correspondent 7 sites d'initiation de la transcription différents !
Une approche possible pourrait être de garder chaque site d'initiation de chaque transcrit pour nos figures, mais cela aboutit à une explosion du nombre de TSS : 200 401 dans la version 27 humaine de GENCODE ! Une autre approche serait de garder uniquement les sites d'initiation uniques, soit environ 7 pour BRCA1. Mais c'est alors faire peser autant un TSS utilisé par des dizaines de modèles de transcrits avec des TSS alternatifs utilisés par de rares transcrits alternatifs. L'approche favorite des bioinformaticiens est donc assez souvent de ne garder qu'un seul site d’initiation de la transcription par gène.
Par facilité, certains s'en sortent en ne gardant que les lignes de type "gene" dans le fichier d'annotation GENCODE (pour plus de détails, voir par exemple cet article). Une ligne par gène, un seul TSS, problème résolu !
Or la ligne gene de GENCODE part du nucléotide le plus en amont parmi tous les transcrit jusqu'au nucléotide le plus en aval parmi tout les transcrits (pas nécessairement dans le même transcrit, d'ailleurs), ce qui nous fera un TSS de BRCA1 ici :
Alors qu'intuitivement, en suivant une règle majoritaire, j'aurais plutôt tendance à privilégier ce site d'initiation là :
Utiliser les lignes "gene" de fichier d'annotation GENCODE pour obtenir une liste d'un unique site d'initiation de la transcription par gène me semble donc assez maladroit pour un nombre important de gènes.
Pour ma part, je prends plutôt le site d'initiation médian de tout les transcrits d'un gène, excluant ainsi les quelques transcrits suspicieux commençant soit très tôt, soit très tard. Dans les cas où il y a un nombre pair de transcrit, je ne prends pas la moyenne des coordonnées des deux TSS du milieu, mais bien l'un des deux (en l’occurrence, le TSS du premier des deux "transcrits du milieu").
Voici par exemple le petit script R que je me suis fait, qui part du fichier GFF3 fournit par GENCODE (pour plus de détails, voir cet article), et retourne un fichier .bed d'un unique TSS par gène. Le script n'est pas des plus rapides malheureusement, mais je n'ai à le faire tourner qu'une seule fois par fichier d'annotation.
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 |
library(rtracklayer) library(tidyverse) # l'import peu prendre 2 minutes en fonction de votre configuration gencode <- import("gencode.v27.annotation.gff3.gz", format = "GFF") %>% # import() retourne un objet GRanges, nous le convertissons en tibble as_tibble transcript <- filter(gencode, type == "transcript") # fonction retournant le transcrit du milieu parmi tous les transcrits du gène 'gene' getMiddleLineFor <- function(gene) { tempt <- dplyr::filter(transcript, gene_id == gene) if(tempt$strand[1] == "+") { # si le gène est '+', le TSS est au 'start' tempt <- dplyr::arrange(tempt, start) } else if(tempt$strand[1] == "-") { # si le gène est '-', le TSS est au 'end' tempt <- dplyr::arrange(tempt, end) } # retourne la ligne du mileur dplyr::slice(tempt, ceiling(nrow(tempt)/2)) } # puis j'applique la fonction à tous les gènes # environ 20 minutes… midTranscript <- map_dfr( unique(transcript$gene_id), getMiddleLineFor ) write.table( select(midTranscript, seqnames, start, end, gene_id,level, strand, gene_type, gene_name), file = "gencode.v27.annotation.middleTSStranscript.bed", quote = FALSE, row.names = FALSE, col.names = FALSE, sep = "\t" ) |
Si je ne veux que la position des TSS et pas le start et le end du transcrit du milieu, je rajoute cette étape :
|
1 2 3 4 5 |
mutate( midTranscript, start = if_else(strand == "+", start, end - 1L), end = if_else(strand == "+", start + 1L, end) ) |
La fonction est facilement adaptable lorsqu'on en souhaite non pas le site d'initiation de la transcription du milieu, mais le site de terminaison de la transcription du milieu !
Après intense réflexion, au lieu de prendre le TSS médian, il pourrait être judicieux de retenir le mode des TSS des transcrits d'un gène ? Il faut cependant une implémentation maline qui traite judicieusement les cas où il y a autant de TSS différents que de transcrits…
Pour aller plus loin
Je vous ai donc proposé une méthode simple pour sélectionner un seul site d'initiation de la transcription par gène. Il est fort probable que pour un certain nombre de gènes, le site d'initiation "majoritaire" soit variable en fonction du tissu : tel TSS sera favorisé dans un type cellulaire, mais un autre TSS le sera dans un autre tissu. Une autre méthode, beaucoup plus lourde, mais plus maline, pourrait donc consister à identifier le TSS majoritaire de vos types cellulaires d’intérêts d'après les données CAGE du consortium Fantom5 (du moins, pour l'humain et la souris).
Une des découvertes de ce consortium a d'ailleurs été que l'idée d'un site d'initiation de la transcription bien défini n'était sans doute pas valide pour tous les gènes. Les promoteurs riches en CpG auraient plutôt une "zone" diffuse d'initiation de la transcription [1].
Enfin, si l'on s’embête tant avec ces histoires de TSS unique par gène, c'est que bien souvent nous ne possédons une mesure d'expression qu'au niveau du gène. Mais les progrès techniques et algorithmiques en RNA-seq font que la quantification au niveau des transcrits directement est de plus en plus fiable [2]. Raisonner avec des transcrits permet d'éviter tous ces détours méthodologiques à la recherche de la "bonne liste" des TSS.
Remerciements
Un grand merci à mes relectrices et relecteurs : Clémence, Nisaea, Yoann M. et Gwenaelle !
Références :
[1] A promoter-level mammalian expression atlas. FANTOM Consortium and the RIKEN PMI and CLST (DGT). DOI : 10.1038/nature13182
[2] Gene-level differential analysis at transcript-level resolution. Lynn Yi, Harold Pimentel, Nicolas L. Bray, Lior Pachter. DOI : 10.1101/190199
Laisser un commentaire