Suivez l'guide :
Manipulation d'intervalles génomiques dans R

Introduction

Nous avons abordé, dans le précédent article de cette série, les bases de la manipulation d'intervalles dans R.

Ce deuxième article a pour objectif de montrer comment manipuler des intervalles génomiques. Au niveau le plus basique, un intervalle est défini par deux nombres entiers positifs délimitant son début et sa fin. Nous allons pouvoir ajouter des couches supplémentaires d'informations, telles que les chromosomes et leurs tailles ainsi que le brin.

Nous allons pour cela utiliser une bibliothèque de la suite Bioconductor, GenomicRanges.

Cette bibliothèque est en quelque sorte l'équivalent des fameux bedtools, très utilisés pour les analyses en génomique.

Installer et charger la bibliothèque GenomicRanges

Commençons par installer et charger la bibliothèque.

Créer un objet GRanges

Coordonnées des intervalles

Nous allons créer un objet contenant des intervalles, ainsi que des informations sur les chromosomes et les brins de chromosomes qui portent ces intervalles. Un objet granges est l'équivalent d'un fichier au format .bed (Browser Extensible Data) contenant au minimum trois colonnes : le nom de la séquence, un intervalle (un objet qui comprend une position de début et une position de fin) et le brin (positif "+" ou négatif "-").

Dans le code ci-dessus, la fonction c() permet de combiner des éléments dans un vecteur (ou une colonne). La fonction rep(x, n) permet de répéter 'n' fois l'objet x. Ainsi, avec le code seqnames = c(rep("chr1", 4), rep("chr2", 4)), nous assignons les quatre premiers intervalles au chromosome 1 (chr1), et les quatre intervalles suivants au chromosome 2 (chr2).

On indique aussi que la moitié de nos intervalles sont sur le brin '+' et le reste sur le brin '-'. Si le brin n'a pas d'importance, on utilise alors strand = "*" lors de la création de l'objet granges. L'argument ignore.strand = TRUE permet d'ignorer le brin dans la plupart des fonctions de la bibliothèque GenomicRanges.

Voyons à quoi ressemble notre nouvel objet :

Métadonnées

Notre objet gr ne contient pour le moment que des informations nécessaire pour la création des intervalles ("GRanges object with 8 ranges and 0 metadata columns"). Nous pouvons ajouter des informations supplémentaires comme des noms de gènes, le taux en GC, etc. à l'aide du paramètre mcols().

Notre objet ressemble maintenant à ceci :

Informations sur les chromosomes

Pour l'instant, les intervalles de notre objet sont situés sur des brins de chromosomes, mais nous n'avons pas précisé leurs tailles, ni le génome ("2 sequences from an unspecified genome; no seqlengths"). Sans ces informations, nous ne pouvons faire des opérations que sur les régions de chromosomes couvertes par ces intervalles.

La commande seqinfo() permet d'afficher et modifier les informations liées aux séquences :

Les informations concernant la taille des chromosomes de l'espèce ou du génome de référence peuvent être trouvées sur le site de l'UCSC ou celui d'Ensembl. La bibliothèque BSgenomes, permet aussi de récupérer ces informations pour certains génomes.

Nous pouvons aussi ajouter ces informations de la façon suivante :

Affichons à nouveau les informations à l'aide de seqinfo() :

Extraire des données

Nous pouvons accéder aux informations à l'aide des commandes de base suivantes :

Nous pouvons également utiliser des opérateurs logiques :

Opérations sur des intervalles génomiques

Un seul objet granges

Nous pouvons effectuer des opérations sur notre objet :

  • reduce(gr) : fusionner les intervalles chevauchants ou adjacents
  • gaps(gr) : extraire les positions non couvertes par des intervalles
  • coverage(gr) : calculer la couverture (nombre de fois que chaque position est couverte par un intervalle)
Résultats des fonctions détaillées ci-dessus

D'autres fonctions existent :

  • flank(gr, width = n, start = TRUE) : extraire une région de longueur 'n' en amont. Pour extraire une région de longueur 'n' en aval, utiliser start = FALSE.
  • shift(gr, shift = n) : déplacer les intervalle de 'n' positions vers la droite (ne tient pas compte du brin). Pour décaler les intervalles vers la gauche, utiliser shift = -n.
  • resize(gr, width = n, fix = "start") : modifier la taille des intervalles de longueur 'n', en déplaçant la fin de l'intervalle (fix = "start") ou le début de l'intervalle (fix = "end").
  • promoters(gr, upstream = x, downstream = y) : extraire les régions promotrices, entre x positions en amont (par défaut 2000) et y positions en aval (par défaut 200)

Plusieurs objets granges

Après avoir vu différentes opérations sur un seul objet granges, voyons comment comparer plusieurs objets.

Commençons par créer deux objets : gr1 et gr2.

Les fonctions suivantes permettent d'extraire l'union et l'intersection de deux objets granges :

  • union(gr1, gr2) : union des deux objets.
  • intersect(gr1, gr2) : intersection des deux objets.

Nous pouvons également rechercher et compter les régions où les intervalles se recouvrent :

  • findOverlaps(gr1, gr2) : détecter les régions où les intervalles se recouvrent. Cette fonction retourne un objet de type Hits qui indique pour chaque intervalle du premier objet (query) les intervalles du deuxième objet qui le recouvrent (subject).
  • countOverlaps(gr1, gr2) pour chaque intervalle du premier objet, indique le nombre d'intervalles du deuxième objet qui le recouvrent.

Pour finir, voyons une dernière fonction qui permet de calculer la distance entre deux objets. Créons tout d'abord deux nouveaux objets :

  • distanceToNearest(gr3, gr4) : calcule la distance entre chaque intervalle du premier objet et les intervalles du deuxième objet.
Résultats des fonctions détaillées ci-dessus

Importer et exporter des fichiers

La bibliothèque rtracklayer permet d'importer et d'exporter des données sous forme d'objet granges.

Attention ! Si vous souhaitez que la colonne de métadonnées soit exportée, elle doit porter l'en-tête name. En effet, les colonnes nommées name et score seront exportées. Pensez bien à vérifier que c'est le cas avant d'exporter votre objet !

Conclusion

Nous avons vu, dans les deux premiers articles de cette série, comment créer et manipuler des intervalles génomiques avec R, à l'aide des bibliothèques IRanges et GenomicRanges.

Dans le prochain article, nous verrons comment effectuer le même genre d'opérations tout en utilisant une syntaxe "tidyverse" à l'aide de la bibliothèque plyranges.

Pour aller plus loin ...


Un grand merci à Guillaume, Matthias et Pierre pour leurs commentaires, conseils et autres remarques constructives !

Les logos des bibliothèques ont été téléchargés sur BiocStickers.

  • À propos de
  • Technicien de la recherche INRAE, travaille sur la diversité génétique du blé tendre, délaissant la paillasse pour l'ordi.

Catégorie: Suivez l'guide | Tags: , , , ,

Laisser un commentaire