Un flux de travail vérifié pour le traitement des données MiRNA-Seq et l’analyse bioinformatique à l’aide de R

Summary

Ici, nous présentons un protocole pour analyser les données de miRNA-Seq à l’aide de R. Le flux de travail permet aux chercheurs d’explorer les réseaux régulés par les miARN et leur importance dans diverses questions biologiques et cliniques. Ce travail vise à servir de guide pratique pour les chercheurs novices et expérimentés dans le domaine de la bioinformatique des miARN.

Abstract

Les microARN (miARN) sont des régulateurs post-transcriptionnels critiques qui influencent un large éventail de processus physiologiques et pathologiques. Avec l’avancement des technologies de séquençage à haut débit, miRNA-Seq est devenu un outil puissant pour profiler les modèles d’expression des miARN. Cependant, l’interprétation fiable de ces données nécessite un pipeline d’analyse standardisé et reproductible. Ici, nous présentons un flux de travail vérifié pour le traitement des données miRNA-Seq et l’analyse bioinformatique à l’aide de R. Ce protocole englobe toutes les étapes essentielles, y compris le prétraitement des données brutes, le contrôle qualité, l’alignement, la quantification, la normalisation, l’analyse d’expression différentielle, la prédiction de cibles, l’enrichissement fonctionnel et la construction de réseaux réglementaires. Conçu pour la flexibilité et la transparence, le flux de travail intègre des packages R largement adoptés et prend en charge l’annotation spécifique à l’espèce et la personnalisation modulaire. De plus, les utilisateurs sont guidés pour effectuer une interprétation biologique en aval en exploitant des bases de données organisées et des outils de visualisation tels que Cytoscape. Ce protocole permet non seulement une analyse statistique robuste, mais aussi des informations significatives sur les interactions miARN-ARNm et leurs rôles dans les mécanismes de la maladie. Il est particulièrement bien adapté aux chercheurs novices et expérimentés qui mènent des études de biomarqueurs de miARN, de modélisation de maladies ou d’études multi-omiques intégratives.

Introduction

Les microARN (miARN) sont de courtes molécules d’ARN non codantes qui influencent significativement l’expression des gènes en agissant au stade post-transcriptionnel¹. Ils fonctionnent généralement en se liant à des séquences complémentaires dans les régions 3' non traduites (UTR) des ARN messagers cibles (ARNm), conduisant à la dégradation de l’ARNm ou à la répression traductionnelle¹. Au cours des deux dernières décennies, les miARN ont été de plus en plus reconnus comme des régulateurs centraux de divers processus biologiques, notamment la prolifération cellulaire, la différenciation, l’apoptose, les réponses immunitaires et le développement des organes². De plus, la dérégulation de l’expression des miARN a été impliquée dans la pathogenèse de nombreuses maladies, telles que le cancer, les maladies cardiovasculaires, les troubles neurologiques et les maladies rénales³. Ces résultats mettent en évidence le potentiel des miARN non seulement en tant que cibles thérapeutiques, mais aussi en tant que biomarqueurs mini-invasifs dans les diagnostics cliniques.

Avec l’avènement des technologies de séquençage de nouvelle génération (NGS), l’étude des miARN est entrée dans une nouvelle ère. Contrairement aux méthodes basées sur des microréseaux qui sont limitées aux miARN connus, le séquençage des miARN (miRNA-Seq) permet un profilage complet, à haut débit et non biaisé des miARN connus et nouveaux dans différents types d’échantillons et conditions⁴. miRNA-Seq offre une sensibilité, une précision et une plage dynamique supérieures, ce qui en fait une méthode privilégiée pour étudier les modèles d’expression des miARN et découvrir les mécanismes de régulation dans des contextes physiologiques et pathologiques⁵. Cependant, l’analyse des données de miRNA-Seq présente des défis de calcul spécifiques, notamment la gestion de courtes longueurs de lecture, la suppression de séquences adaptatrices, la distinction entre les membres de la famille des miARN étroitement apparentés et la gestion d’une redondance élevée dans le nombre de lectures⁶. Ces caractéristiques nécessitent un flux de travail analytique soigneusement conçu et standardisé.

Bien que divers pipelines et outils logiciels aient été développés pour l’analyse des données miRNA-Seq, beaucoup d’entre eux reposent sur des interfaces utilisateur graphiques ou des flux de travail fixes qui limitent la flexibilité et la reproductibilité⁷. En revanche, l’environnement de programmation R fournit une plate-forme puissante et personnalisable pour l’analyse bio-informatique⁸. R offre un riche écosystème de packages pour la modélisation statistique, la visualisation des données et l’intégration avec des bases de données biologiques. Cela permet aux utilisateurs d’effectuer des analyses complètes et reproductibles de manière transparente et basée sur des scripts. De plus, la nature modulaire des flux de travail R permet aux chercheurs d’adapter chaque étape en fonction d’exigences expérimentales spécifiques, du prétraitement des données brutes à l’interprétation fonctionnelle.

Dans ce protocole, nous présentons un flux de travail d’analyse miRNA-Seq vérifié et complet entièrement mis en œuvre dans R, dans le but de fournir une solution reproductible et adaptable à l’utilisateur pour les chercheurs travaillant avec des données d’expression de miARN. Le flux de travail commence par le contrôle de la qualité et le découpage de l’adaptateur des lectures de séquençage brutes, suivis de l’alignement sur un génome de référence ou des séquences de miARN connues. Les étapes suivantes comprennent la quantification du nombre de lectures, la normalisation, l’analyse de l’expression différentielle, la prédiction des gènes cibles, l’enrichissement fonctionnel et la visualisation du réseau. Le flux de travail intègre plusieurs packages R largement utilisés et bien entretenus, garantissant à la fois la fiabilité et la compatibilité avec les futures mises à jour et extensions.

L’une des principales forces de ce protocole réside dans sa capacité à aller au-delà des résultats d’expression différentielle et à fournir une interprétation biologique significative. En intégrant des bases de données organisées d’interactions miARN-ARNm validées et prédites, le flux de travail permet aux utilisateurs d’identifier des gènes cibles biologiquement pertinents. Ces cibles peuvent ensuite être soumises à des analyses d’ontologie génétique et d’enrichissement des voies pour découvrir les processus biologiques et les voies moléculaires affectés. Dans la dernière étape, les réseaux d’interaction miARN-ARNm peuvent être visualisés à l’aide d’outils externes tels que Cytoscape⁹, ce qui donne un aperçu du paysage réglementaire et identifie les principaux miARN pivots ayant une importance fonctionnelle potentielle.

Cette méthode a été appliquée avec succès dans des contextes de recherche clinique, y compris des études sur les maladies rénales, où les miARN circulants servent de biomarqueurs prometteurs pour le diagnostic et le pronostic¹⁰. Cependant, la conception modulaire et flexible du flux de travail le rend adapté à un large éventail d’applications, notamment la modélisation de maladies, les études de réponse aux médicaments, la biologie du développement et la génomique comparative. Les chercheurs peuvent facilement adapter le flux de travail pour prendre en charge des annotations spécifiques à une espèce, des conditions expérimentales ou des couches supplémentaires de données omiques.

En offrant une solution open-source basée sur des scripts, ce pipeline centré sur R répond à plusieurs des principales limitations associées aux outils miRNA-Seq existants, notamment la personnalisation limitée, la dépendance à des interfaces graphiques non transparentes, le manque de prise en charge des organismes non modèles, la faible reproductibilité due à l’absence de contrôle de version et la difficulté à s’intégrer aux cadres d’analyse statistique et fonctionnelle en aval. Il permet un contrôle total des paramètres de traitement des données, encourage la reproductibilité grâce à un code contrôlé par version et favorise la transparence dans la recherche bioinformatique. Alors que l’importance des miARN ne cesse de croître dans le contexte de la biologie des systèmes et de la médecine translationnelle, il devient de plus en plus essentiel d’avoir accès à un cadre d’analyse fiable et adaptable.

Protocol

REMARQUE : Les matériaux avec des liens vers des logiciels sont répertoriés dans la table des matériaux.

1. Préparer des échantillons d’ARN et des banques de séquences

REMARQUE : Effectuez l’extraction et le séquençage de l’ARN en dehors de ce flux de travail informatique. Il existe plus d’une façon d’analyser les données de séquençage des miARN. Cette section fournit le contexte d’une pratique.

1. Extraire l’ARN total : Extraire l’ARN total des échantillons biologiques à l’aide d’un kit optimisé pour l’isolement de petits ARN (par exemple, un kit d’isolement de miARN). Suivez attentivement le protocole du fabricant. Assurez-vous d’utiliser des consommables sans RNase et conservez les échantillons sur de la glace pour minimiser la dégradation.
2. Évaluer l’intégrité et la quantité d’ARN : Exécutez 1 à 2 μL de l’ARN extrait sur un bioanalyseur ou un appareil équivalent. Vérifiez le numéro d’intégrité de l’ARN (RIN) et assurez-vous qu’il est ≥ 7.0 pour un séquençage fiable. Enregistrez la concentration à l’aide d’un spectrophotomètre ou d’un fluorimètre.
3. Construire de petites banques d’ARN : Préparez des banques de séquençage à partir de 1 μg d’ARN total à l’aide d’un kit de préparation de petites banques d’ARN-seq. Suivez le protocole du kit pour ligaturer les adaptateurs, transcrire inverser et amplifier l’ADNc. Purifiez les produits PCR avec une sélection de taille (par exemple, insert de 18 à 30 nt) pour enrichir les fragments de miARN.
4. Bibliothèques de séquences : chargez les bibliothèques sur une plate-forme de séquençage à haut débit. Définissez la configuration d’exécution pour le séquençage à une extrémité avec une longueur de lecture de ~50 pb. Assurez-vous que chaque échantillon génère environ 10 millions de lectures brutes pour obtenir une profondeur suffisante.
5. Exporter la sortie du séquençage : Après le séquençage, exportez les données brutes sous forme de fichiers FASTQ à l’aide du logiciel de sortie de données de l’instrument. Vérifiez que le répertoire de sortie contient à la fois les lectures de séquence et les fichiers de score de qualité associés. Stockez les fichiers FASTQ dans un répertoire structuré pour l’analyse en aval.

2. Prétraitez les lectures brutes et effectuez le contrôle qualité

1. Séquences d’adaptateur de trim
   1. Installez et configurez Cutadapt ou fastp.
   2. Exécutez le découpage de l’adaptateur sur chaque fichier FASTQ à l’aide de la commande suivante :
      cutadapt -a XXXX -o trimmed_reads.fastq raw_reads.fastq
      REMARQUE : La séquence d’adaptateur suivie de « -a » doit s’appliquer avec le kit de préparation de banque d’ARNs spécifique. '-o' définit le nom du fichier de sortie suivi du fichier d’entrée.
2. Évaluer la qualité de lecture
   1. Utilisez FastQC pour générer des rapports de contrôle de la qualité :
      fastqc trimmed_reads.fastq
   2. Examiner les scores de qualité par base, la distribution de la longueur de lecture et la contamination de l’adaptateur : Ouvrez le rapport HTML FastQC généré pour chaque fichier FASTQ dans un navigateur Web. Examinez les modules suivants étape par étape :
      1. Qualité de la séquence par base : Assurez-vous que la plupart des bases se trouvent dans la zone verte (score Phred ≥30). Recherchez toute baisse de qualité à l’extrémité 3' qui pourrait indiquer des erreurs de séquençage.
      2. Distribution de la longueur de lecture : Confirmez que la distribution correspond à la taille d’insert attendue (par exemple, 18-30 nt pour les miARN). Vérifiez qu’il n’y a pas de pics inattendus.
      3. Contenu de l’adaptateur : vérifiez que les séquences d’adaptateur ont été correctement tronquées. Vérifiez que le pourcentage de contamination de la carte est proche de zéro après le rognage.
   3. Enregistrez le rapport récapitulatif FastQC et signalez tous les échantillons dont les mesures sont de mauvaise qualité pour les réajuster ou les exclure d’une analyse plus approfondie.

3. Cartographiez les lectures et générez des matrices de comptage

1. Aligner les lectures sur la référence
   1. Téléchargez le fichier FASTA du génome de référence ou des séquences de miARN matures (par exemple, à partir de miRBase)¹¹. Exemple:
      wget ftp://mirbase.org/pub/mirbase/CURRENT/mature.fa
   2. Indexez le génome de référence à l’aide de Bowtie.
      1. Ouvrez un terminal et exécutez la commande suivante pour générer l’index :
         bowtie-build reference.fa reference_index
      2. Remplacez reference.fa par le nom de fichier FASTA réel.
      3. Remplacez reference_index par le préfixe souhaité pour l’index.
      4. Assurez-vous que Bowtie génère plusieurs fichiers d’index (par exemple, .ebwt). Vérifiez que ces fichiers sont présents dans le répertoire de travail, car ils sont nécessaires à l’alignement.
   3. Alignez les lectures à l’aide du nœud papillon avec les paramètres appropriés pour les lectures courtes. Exemple:
      bowtie -v 0 -a --best --strata reference_index trimmed_reads.fastq > aligned_reads.sam
      REMARQUE : Le fichier d’entrée est trimmed_reads.fastq et le fichier de sortie est aligned_reads.sam. '-v 0' signifie que nous ne permettons aucune erreur dans toute la lecture. '-a -best -strata' signifie que nous jetons tout alignement qui a plus de décalages que le meilleur.
2. Quantifier l’expression des miARN
   1. Convertissez les fichiers SAM au format BAM à l’aide de SAMtools.
      samtools view -S -b aligned_reads.sam > aligned_reads.bam
      REMARQUE : Le fichier d’entrée aligned_reads.sam est le résultat de la dernière commande. Le fichier de sortie aligned_reads.bam est préparé pour l’analyse ultérieure.
      1. Utilisez SAMtools pour compresser et trier le fichier d’alignement :
         samtools sort aligned_reads.bam -o aligned_reads_sorted.bam
samtools index aligned_reads_sorted.bam
      2. Assurez-vous que la première commande convertit le fichier SAM au format BAM.
      3. Assurez-vous que la deuxième commande trie le fichier BAM par coordonnées génomiques.
      4. Assurez-vous que la troisième commande génère un fichier d’index (.bai), qui est requis pour les analyses en aval.
      5. Vérifiez que le fichier BAM trié et son index ont été créés avec succès avant de passer à la quantification.
   2. Utilisez featureCounts ou HTSeq-count pour générer une matrice de comptage à l’aide de l’annotation GTF d’annotation de miARN :
      featureCounts -a miRNA.gtf -o counts.txt aligned_reads.bam
      REMARQUE : featureCounts quantifie les lectures dans le fichier d’entrée aligned_reads.bam basé sur miRNA.gtf, et les sorties counts.txt.

4. Effectuer une analyse d’expression différentielle dans R

1. Données de comptage de charges
   1. Importez les métadonnées de la matrice de comptage et de l’échantillon dans R :
      library(DESeq2)
countData <- read.csv("counts.csv", row.names=1)
colData <- read.csv("metadata.csv", row.names=1)
dds <- DESeqDataSetFromMatrix(countData = countData, colData = colData, design = ~ condition)
      REMARQUE : Les données de comptage (counts.csv) et le groupe d’échantillons (metadata.csv) doivent être fournis à DESeq2. 'condition' clarifie ici le groupe d’échantillons fournis. Pour les besoins individuels, reportez-vous au manuel de DESeq2¹².
2. Normaliser et transformer les données
   1. Normalisez les données de comptage à l’aide de la méthode par défaut de DESeq2 :
      dds <- DESeq(dds)
   2. Effectuer une transformation de stabilisation de la variance :
      vsd <- vst(dds, blind=FALSE)
   3. Visualisez le clustering d’échantillons à l’aide de l’ACP :
      plotPCA(vsd, intgroup="condition")
3. Identifier les miARN exprimés de manière différentielle
   1. Extraire et trier les résultats de l’expression différentielle :
      res <- results(dds)
resOrdered <- res[order(res$pvalue), ]
      REMARQUE : Nous réorganisons le fichier de résultats 'res' en fonction de la valeur de pvalue.
      summary(res)
   2. Filtre pour les miARN exprimés de manière significativement différentielle (valeur p < 0,05, |log2FC| > 1).
      sig_miRNA <- subset(res, pvalue < 0.05 & abs(log2FC) > 1)
      REMARQUE : Il existe plusieurs seuils pour filtrer les miARN considérablement modifiés. Les exigences 'valeur p < 0,05, |log2FC| > 1' sont largement appliqués. Les seuils peuvent être ajustés pour des données individuelles.
4. Visualiser les modifications d’expression
   1. Installez et chargez le paquet EnhancedVolcano.
   2. Créez le graphique du volcan :
      library(EnhancedVolcano)
EnhancedVolcano(res,
lab = rownames(res),
x = 'log2FoldChange',
y = 'pvalue',
title = 'Differentially Expressed miRNAs')

5. Prédire les gènes cibles des miARN

1. Interroger des bases de données
   1. Utilisez des ressources en ligne telles que TargetScan, miRDB et miRTarBase pour rechercher des microARN spécifiques et récupérer des gènes cibles.
   2. Concentrez-vous sur des cibles validées expérimentalement pour une plus grande confiance.
2. Automatiser la prédiction dans R
   1. Chargez le package multiMiR et interrogez les cibles validées :
      library(multiMiR)
target_results <- get_multimir(mirna = c("hsa-miR-21-5p"), table = "validated")
      REMARQUE : Ici, nous prenons « has-miR-21-5p » par exemple et obtenons des cibles validées.
   2. Extraire des symboles de gènes cibles uniques pour l’analyse de l’enrichissement :
      genes <- unique(target_results@data$target_symbol)

6. Effectuer une analyse d’enrichissement fonctionnel

1. Effectuer l’enrichissement GO
   1. Outils d’enrichissement de charge :
      library(clusterProfiler)
library(org.Hs.eg.db)
      REMARQUE : Ici, nous chargeons une base de données contenant des annotations du génome humain, qui sont utiles pour convertir les identificateurs de gènes courants.
   2. Exécutez l’analyse d’enrichissement GO pour les procédés biologiques :
      ego <- enrichGO(gene = genes,
OrgDb = org.Hs.eg.db,
keyType = "SYMBOL",
ont = "BP",
pAdjustMethod = "BH",
pvalueCutoff = 0.05)
dotplot(ego)
      REMARQUE : Pour l’utilisation d’enrichGO, nous fournissons une liste de gènes, et précisons qu’ils sont nommés comme « SYMBOLE » ici. Nous avons réalisé un enrichissement par procédé biologique, qui correspond au paramètre 'ont ='BP'. Nous effectuons plusieurs corrections de tests, nous désignons donc pAdjustMethod = « BH ». Pour un seuil significatif, nous choisissons pvalueCutoff = 0.05.dotplot présente un résultat visualisé. Pour plus d’options individuelles, reportez-vous au manuel de clusterProfiler¹³.
2. Effectuer l’enrichissement de la voie KEGG
   1. Exécuter l’enrichissement KEGG :
      ekegg <- enrichKEGG(gene = genes, organism = 'hsa')
dotplot(ekegg)
      REMARQUE : Pour l’utilisation d’enrichKEGG, nous fournissons une liste de gènes et précisons que l’organisme est humain (« a »). dotplot présente un résultat visualisé. Pour plus d’options individuelles, reportez-vous au manuel de clusterProfiler¹³.

7. Construire et visualiser le réseau d’interaction miARN - ARNm

1. Exporter des données pour la visualisation du réseau
2. Créez une trame de données de paires de gènes miARN-cible basée sur les cibles générées à partir de TargetScan, miRDB ou miRTarBase.
3. Écrivez la table réseau au format CSV :
   write.csv(miRNA_target_pairs, "network.csv")
4. Importer dans Cytoscape
   1. Ouvrez Cytoscape et importez la table du réseau.
   2. Visualisez le réseau à l’aide d’une disposition dirigée par la force ou circulaire.
   3. Analysez les propriétés topologiques (par exemple, le degré de centralité) pour identifier les miARN pivots.

Results

Nous avons téléchargé la matrice d’expression des microARN à partir de GSE133530 et effectué directement une analyse d’expression différentielle. Nous avons fourni un exemple de script R analytique pour l’ensemble de données dans le fichier supplémentaire 1. L’ensemble de données a effectué un profilage global des miARN sur 16 kystes rénaux de différentes tailles (kystes minimes : moins de 1 à 5 ml, n = 10 ; kystes moyens : entre 10 et 25 ml, n = 4 ; gros kystes : supérie...

Discussion

L’analyse des données miRNA-Seq présente des défis distincts en raison de la petite taille et de la redondance des lectures, ce qui rend essentiel un contrôle qualité rigoureux et un prétraitement. L’une des étapes les plus importantes du flux de travail est le découpage de l’adaptateur. Étant donné que les miARN ont une longueur d’environ 22 nucléotides, les séquences adaptatrices peuvent facilement dominer les lectures si elles ne sont pas correctement éliminées. L’absence d’un réglage...

Materials

List of materials used in this article

Name	Company	Catalog Number	Comments
Agilent-021827 Microréseau de miARN humain	Agilent	/	Une puce commerciale pour le profilage des microARN d’échantillons humains
Nœud papillon	Université Johns Hopkins	http://bowtie-bio.sourceforge.net/index.shtml	Un outil logiciel pour aligner les lectures de séquençage sur les longues séquences de référence
clusterProfiler (package R)	Bioconducteur	https://bioconductor.org/packages/clusterProfiler/	Un package R conçu pour l’analyse d’enrichissement fonctionnel et la visualisation de données biologiques à haut débit.
Cutadapt	Libre	https://cutadapt.readthedocs.io	Outil en ligne de commande qui supprime les séquences d’adaptateur, les amorces, les queues poly-A et autres fragments indésirables des lectures de séquençage à haut débit.
Cytoscape	Cytoscape Consortium	https://cytoscape.org/	Une plate-forme logicielle open-source conçue pour la visualisation et l’analyse de réseaux biologiques complexes.
DESeq2 (package R)	Bioconducteur	https://bioconductor.org/packages/DESeq2/	Un package R conçu pour l’analyse différentielle de l’expression génique des données de comptage
EnhancedVolcano (forfait R)	Bioconducteur	https://bioconductor.org/packages/EnhancedVolcano/	; Un package R conçu pour créer des graphiques de volcans de qualité publication.
FastQC	Babraham Bioinformatique	https://www.bioinformatics.babraham.ac.uk/projects/fastqc/	Un outil de contrôle qualité open source pour les données de séquençage à haut débit.
featureCounts	Sous-lecture / SourceForge	http://subread.sourceforge.net/	Programme utilisé pour compter les lectures associées à des caractéristiques génomiques
Nombre de HTSeq	Paquet Python	https://htseq.readthedocs.io	Outil en ligne de commande qui compte le nombre de lectures de séquençage à haut débit alignées qui chevauchent des caractéristiques génomiques telles que des gènes ou des exons. Je
Puce d’expression de microARN Illumina Human v2	Illumina  ;	/	Une puce commerciale pour le profilage des microARN d’échantillons humains
multiMiR (package R)	Bioconducteur	https://bioconductor.org/packages/multiMiR/	Un package R qui fournit la plus grande collection intégrée de microARN prédits et validés expérimentalement. Cibler les interactions ainsi que leurs associations avec des maladies et des médicaments.
Org. Hs.eg.db (package R)	Bioconducteur	https://bioconductor.org/packages/org.Hs.eg.db/	Un progiciel d’annotation conçu pour la recherche en génomique humaine (Homo sapiens).
Logiciel R	Projet R	https://www.r-project.org/	Un projet open-source pour le calcul statistique
Rstudio	Posit PBC	/	Un environnement de développement intégré permet d’être plus productif avec R et Python
SAMtools	Libre	http://www.htslib.org/	Progiciel de manipulation des données de séquençage de nouvelle génération (NGS).