Avec l’avènement des nouvelles technologies d’analyse du transcriptome, comme les puces à ADN, il est devenu de plus en plus aisé d’obtenir un portrait de l’expression de tous les gènes d’un individu à un temps donné. Les plateformes de chips, comme celles de Affymétrix ou de Illumina (les deux compagnies qui ont la plus grosse part du marché), permettent d’estimer le niveau d’expression de tous les gènes connus chez l’homme en une seule journée de manipulations.
Ces expériences génèrent toutefois un nombre impressionnant de données et il devient parfois difficile de dégager les conclusions biologiques conséquentes avec le design expérimental. Dans ce projet de bio-informatique, nous avons traité un jeu de données nommé GSE61276 afin d’y réaliser : des contrôles de qualité, des analyses différentielles et une analyse des voies métaboliques les plus représentées.
Le jeu de données a été généré avec une puce GPL10558 “Human HT-12 v4” Expression BeadChip de Illumina. Ce modèle de chips compte environ 47 000 sondes ciblant environ 31 000 transcrits et transcrits alternatifs de gènes connus chez l’homme. L’équipe de recherche biomédicale qui a généré les données ont publié leur conclusions dans la revue BMC Genomics en 2014 sous le titre : “Genetic and epigenetic regulation of gene expression in fetal and adult human livers”. par : Marc Jan Bonder etal. (Estonian Genome Center, University of Tartu).
Pour leur analyse, ils ont prélevé 106 échantillons de biopsies de foie provenant de trois groupes: fœtal (14 échantillons), donneurs d’organes (50 échantillons) et ablation partielle du foie (42 échantillons). Dans leur article, les chercheurs ont aussi fait des analyses de méthylation (qui ne seront pas traités ici) et ont entrecroisé leurs résultats d’expression et de méthylation. Ils ont trouvé que ce qui module le plus ces profils d’expression est la différence d’âge et le stade de développement des tissus hépatiques. Ils ont identifié 657 gènes méthylés différentiellement chez l’adulte, s’apparentant aux facteurs de transcription HNF1A (Hepatocyte Nuclear Factor 1 Alpha) et HNF4A (Hepatocyte Nuclear Factor 4 Alpha) et environ 1 000 gènes spécifiques au foie fœtal, enrichis pour les sites de liaison GATA1 (Erythroid transcription factor), STAT5A (Signal Transducer and Activator of Transcription 5A), STAT5B (Signal Transducer and Activator of Transcription 5B) et YY1 (Yin Yang 1). Ces conclusions sont peu surprenantes, alors que les HNF sont des transcrits associés à la différenciation des cellules hépatiques. Les gènes STAT (facteur de transcription) et YY1 (méthylation) sont impliqués dans la régulation de la différenciation cellulaire. Pour ce qui est de GATA1, il est possible que le fait de l’observer dans des biopsies de tissus hépatique a un lien avec le rôle du foie dans l’élimination des cellules sanguines (leucocytes et mégalocytes) usées.
Dans notre traitement des données, une analyse différentielle a aussi été réalisée entre chaque paire de groupes présents. En plus d’analyser les patrons d’expression entre les cellules hépatiques fœtales et adultes, le but de ces analyses différentielles est de voir s’il est possible de dégager des similitudes entre les échantillons provenant de d’ablation de lobes de foie atteints d’hépatocarcinomes (les portions non-tumorales) avec ceux de biopsies provenant de dons d’organes (des foies présumés sains).
Les analyses statistiques ont été réalisées avec RStudio et différents packages d’analyse de data de micro-array : “GOstats” , “hgu133plus2.db” , “limma” , “lumi” , “lumiHumanAll.db” ,“lumiHumanIDMapping” , “biomaRt” , “pcaMethods” ,“KEGGprofile” , “knitr” et “GOstats”.
Normalisation des données
Les données ont été normalisées avec la fonction RSN (pour Robust Spline Normalization) afin de rendre les données comparables entre-elles. En effet, un boxplot des données non-normalisées montrait clairement que certaines données n’étaient pas comparables aux autres. Cette normalisation a donné des résultats qui semblaient bons, du moins selon les boxplots et les PCA.
Boxplot et PCA sur les données normalisées par RSN. Le groupe en bleu représente les échantillons fœtaux, en noir Adult_Accident et en rouge Adult_Cancer
Toutefois, comme le montre les MA-plots et les volcano-plots, les données ont été normalisées d’une drôle de façon, avec des données ayant un FC beaucoup trop élevé. Il aurait possible de spécifier un seuil de « cut-off » lors de la normalisation (par exemple avec « excludeFold = 2 »), mais cela aurait entrainé la perte de signal biologique intéressant.
Les auteurs de l’étude (Bonder et al. 2014) ont aussi mis en ligne un jeu de données normalisées où les boxplot sont tous à des hauteurs comparables. Il semble que ce jeu de données ne comprenne que environ 21 000 gènes, contre plus de 47 000 gènes pour le jeu de données complet, mais les résultats biologiques sont plus intéressants pour la suite du pipeline.
Analyses différentielles
Afin de conduire des analyses différentielles et de pouvoir interpréter les résultats obtenus plus facilement, la matrice correspondant aux gènes a été annotée selon la base de données de Illumina (dbSpecie <- “lumiHumanAll.db”) selon le génome de l’homo sapiens (specie <- “hsa”) et le système d’ID de Illumina (IDMapping <- “lumiHumanIDMapping”). Aussi directement sur la matrice « lumi », les nuID ont été ajoutés [ addNuID2lumi(lumi, lib.mapping= IDMapping) ]. Pour être exact, ces annotations ont été greffées à la matrice d’expression dès son importation et avant les tentatives de normalisation.
Une fois ces annotations faites, d’autres informations sont obtenues grâce à biomaRt. Les attributs “hgnc_symbol”, “entrezgene”, “ensembl_transcript_id” et “description” sont obtenus pour tous les gènes présents sur la chip. Cela permet de relier les ID des sondes aux gènes correspondant dans la base de données de Ensembl, qui possède un système très complet d’annotations de gènes.
Les analyses différentielles sont ensuite conduites avec la fonction lmFit et la fonction eBayes du package Lumi. Les comparaisons sont ensuite réalisées pour chaque paire de conditions possibles (Fetal, Adulte_Cancer, Adulte_Accident). Les résultats sont aussi présentés sous la forme de diagrammes de Venne, où il est possible de voir facilement le nombre de gènes communs aux paires de groupes. Des histogrammes simples des valeurs de p (axe des x = p-value et axe des y = nombre d’échantillons dans cet intervalle).
Le package KEGG a ensuite été utilisé pour faire une analyse des voies métaboliques. Avec les tables des gènes les plus modulés en entrée, les fonctions du package ont permis de générer des tableaux des Ontologies de Gènes les plus représentées. De plus, il a été possible de générer des cartes des voies métaboliques les plus représentées, détaillées au niveau de tous les métabolites impliqués (protéines, molécules activatrices, substrats, etc.)
Les résultats sont présentés dans le pipeline ProjetStats.rmd:
http://rpubs.com/jauger/79243
Dans un premier temps, nous avons normalisé les données brutes disponibles, mais il n’a pas été possible d’obtenir un résultat satisfaisant pour la suite des analyses. Le reste de l’analyse a été réalisé sur un jeu de données préalablement normalisé par les auteurs de l’étude.
Selon la PCA, les deux premières composantes expliquent environ la moitié de la variation des échantillons ensemble. Les grappes distinctes formées par les données montrent que les groupes semblent avoir des profils d’expression qui leur sont propres. La majorité de la variance semble expliquée par la composante 1, ce qui suggère qu’il y a une explication de la variance par les conditions expérimentales. La séparation des couleurs selon les groupes suggère que les groupes expérimentaux ont un impact sur les profils d’expression.
Les M versus A plots montrent distribution assez normale. Il y a quelques points qui s’écartent du centre, mais il n’y a pas de tendance assez marquée pour indiquer un défaut des données utilisées.
Les résultats affichés dans les topTables annotées sont difficiles à interpréter en soi. Il est intéressant d’analyser les résultats d’ontologies, qui els résume à leur façon.
Des six voies mises de l’avant par ce tableau, il y en a trois qui ont rapport aux virus et trois autres au transport cellulaire et interactions entre cellules d’un organisme multicellulaire. Ces voies sont tout à fait conséquentes avec ce à quoi on pourrait s’attendre, alors que la principale différence entre le groupe fœtal et le groupe adulte est le fait que le foie soit en développement ou entièrement actif. Les gènes associés aux virus sont plus fortement actifs chez l’adulte, alors que son organisme est constamment en train de se défendre contre les pathogènes, alors que le fœtus est dans un environnement entièrement aseptisé. Les trois autres gènes modulés sont des gènes qui sont associés à la communication intercellulaire et au développement. Il n’est pas surprenant que ces gènes soient modulés dans un organisme en plein développement.
Ensuite, la deuxième paire de groupes comparés au niveau de l’ontologie est constituée des Adultes-Cancer V.S. Adultes-Accident (présumés sains au niveau du cancer). Quatre voies ont rapport à la manipulation de l’ADN nucléaire ou des transcrits d’ARN. Ces voies métaboliques sont donc différentiellement exprimées entre des échantillons de tissus sains et des tissus associés à des cancers. Ces gènes sont donc des candidats potentiels pour la compréhension des mécanismes d’oncogenèse dans les tissus hépatiques.
Le schéma des voies métaboliques proposé par KEGG met de l’avant la voie appelée « voie de la caféine ». Cela s’explique bien par rapport aux conditions expérimentales. En effet, il faut considérer que les échantillons comparés sont des échantillons Adultes V.S. Fœtaux. Comme mentionné plus tôt, les foies ont des profils d’expression qui leur sont très spécifiques, alors que ces profils sont modulés à la fois par le stade de développement (l’âge du patient) et par une variation individuelle (qui pourrait être expliquée par l’environnement biotique et abiotique du patient). En effet, il n’est pas surprenant que les profils d’expression des foies fœtaux ne soient pas en réaction à des agents xénobiotiques, comme la caféine. De plus, les protéines surlignées par cette voie métaboliques sont des Cytochromes (CYP1A2, CYP2A6, etc.) qui sont donc impliqués dans la voie métabolique plus large de la détoxification hépatique de l’organisme.
En fait le fait que ce soit le métabolisme de la caféine qui ressorte en premier est signe du fait que la différence entre les groupes est surtout expliquée par le stade de développement (et donc de l’exposition aux xénobiotiques). Cela veut aussi dire qu’il est difficile de différencier les autres signaux biologiques du bruit engendré par les composés en cours de détoxification. Il serait intéressant d’obtenir plus d’informations par rapport aux habitudes alimentaires et de consommation des participants afin de les relier aux profils d’expression étudiés. Par exemple, s’il était possible d’associer un signal biologique aux habitudes de consommation de café, il serait éventuellement possible de faire un modèle linéaire robuste. Ce modèle pourrait par la suite permettre d’enlever le signal dû à la caféine et ainsi pouvoir observer d’autres voies métaboliques en action dans le foie.
La seconde voie métabolique proposée par KEGG est celle de la cascade complément-coagulation. Cette conclusion va un peu de pair avec la conclusion de l’équipe de Bonder comme quoi les protéines comme GATA1 (Erythroid transcription factor) sont surexprimées dans le foie adulte. En effet,, comme mentionné plus tôt, le foie joue un rôle primordial dans l’élimination des cellules sanguines usées. Il est très possible que des cellules sanguines se soient retrouvées dans les biopsies t analysées avec les cellules hépatiques. Il n’y a pas de mentions de précautions à cet égard dans le protocole de Bonder. La voie du complément (cascade enzymatique visant la destruction de parasites) et la coagulation sanguine sont très clairement des voies métaboliques associées aux fonctions des cellules sanguines.
En somme, il est possible de dégager des conclusions biologiques des données présentées dans le jeu de données GSE61276 et de les expliquer. Toutefois, ces conclusions sont surtout relatives au stade développemental et à la variation inter-sujets. Il serait intéressant de tenter d’enlever ces “bruits biologiques” pour pouvoir étudier les pathologies comme le cancer. Il serait aussi intéressant d’obtenir des résultats de données d’expression afin de les recouper avec nos données.
[1] B. M. Bolstad, B. M. Bolstad, R. a Irizarry, R. a Irizarry, M. Strand, M. Strand, T. P. Speed, and T. P. Speed, “A comparison of normalization metholds for high density oligonucleotide array data based on variance and bias,” Bioinformatics, vol. 19, no. 2, pp. 185–193, 2003.
[2] G. Courtois, J. G. Morgan, L. A. Campbell, G. Fourel, and G. R. Crabtree, “Interaction of a liver-specific nuclear factor with the fibrinogen and alpha 1-antitrypsin promoters.,” Science, vol. 238, no. 4827, pp. 688–92, Oct. 1987.
[3] A. Gerrits, Y. Li, B. M. Tesson, L. V Bystrykh, E. Weersing, A. Ausema, B. Dontje, X. Wang, R. Breitling, R. C. Jansen, and G. de Haan, “Expression quantitative trait loci are highly sensitive to cellular differentiation state.,” PLoS Genet., vol. 5, no. 10, p. e1000692, Oct. 2009.
[4] W. Li, “Volcano plots in analyzing differential expressions with mRNA microarrays.,” J. Bioinform. Comput. Biol., vol. 10, no. 6, p. 1231003, Dec. 2012.
[5] M. J. Bonder, S. Kasela, M. Kals, R. Tamm, K. Lokk, I. Barragan, W. A. Buurman, P. Deelen, J.-W. Greve, M. Ivanov, S. S. Rensen, J. V van Vliet-Ostaptchouk, M. G. Wolfs, J. Fu, M. H. Hofker, C. Wijmenga, A. Zhernakova, M. Ingelman-Sundberg, L. Franke, and L. Milani, “Genetic and epigenetic regulation of gene expression in fetal and adult human livers.,” BMC Genomics, vol. 15, p. 860, Jan. 2014.