Un flusso di lavoro verificato per l'elaborazione dei dati MiRNA-Seq e l'analisi bioinformatica utilizzando R

Summary

Qui, presentiamo un protocollo per analizzare i dati miRNA-Seq utilizzando R. Il flusso di lavoro consente ai ricercatori di esplorare le reti regolate dai miRNA e il loro significato in diverse questioni biologiche e cliniche. Questo lavoro intende fungere da guida pratica sia per i ricercatori alle prime armi che per quelli esperti nel campo della bioinformatica dei miRNA.

Abstract

I microRNA (miRNA) sono regolatori post-trascrizionali critici che influenzano un'ampia gamma di processi fisiologici e patologici. Con il progresso delle tecnologie di sequenziamento ad alto rendimento, il miRNA-Seq è emerso come un potente strumento per la profilazione dei modelli di espressione dei miRNA. Tuttavia, un'interpretazione affidabile di tali dati richiede una pipeline di analisi standardizzata e riproducibile. Qui, presentiamo un flusso di lavoro verificato per l'elaborazione dei dati miRNA-Seq e l'analisi bioinformatica utilizzando R. Questo protocollo comprende tutte le fasi essenziali, tra cui la pre-elaborazione dei dati grezzi, il controllo qualità, l'allineamento, la quantificazione, la normalizzazione, l'analisi dell'espressione differenziale, la previsione del target, l'arricchimento funzionale e la costruzione della rete regolatoria. Progettato per garantire flessibilità e trasparenza, il flusso di lavoro integra pacchetti R ampiamente adottati e supporta l'annotazione specifica per specie e la personalizzazione modulare. Inoltre, gli utenti sono guidati a condurre l'interpretazione biologica a valle sfruttando database curati e strumenti di visualizzazione come Cytoscape. Questo protocollo non solo supporta una solida analisi statistica, ma consente anche di ottenere informazioni significative sulle interazioni miRNA-mRNA e sui loro ruoli nei meccanismi della malattia. È particolarmente adatto sia per i ricercatori alle prime armi che per quelli esperti che conducono la scoperta di biomarcatori di miRNA, la modellazione di malattie o studi multi-omici integrativi.

Introduction

I microRNA (miRNA) sono molecole di RNA corte non codificanti che influenzano significativamente l'espressione genica agendo nella fase post-trascrizionale¹. Tipicamente funzionano legandosi a sequenze complementari nelle regioni 3' non tradotte (UTR) degli RNA messaggeri bersaglio (mRNA), portando alla degradazione dell'mRNA o alla repressione traduzionale¹. Negli ultimi due decenni, i miRNA sono stati sempre più riconosciuti come regolatori centrali di vari processi biologici, tra cui la proliferazione cellulare, la differenziazione, l'apoptosi, le risposte immunitarie e lo sviluppo degli organi². Inoltre, la disregolazione dell'espressione dei miRNA è stata implicata nella patogenesi di numerose malattie, come il cancro, le malattie cardiovascolari, i disturbi neurologici e le malattie renali³. Questi risultati evidenziano il potenziale dei miRNA non solo come bersagli terapeutici, ma anche come biomarcatori minimamente invasivi nella diagnostica clinica.

Con l'avvento delle tecnologie di sequenziamento di nuova generazione (NGS), lo studio dei miRNA è entrato in una nuova era. A differenza dei metodi basati su microarray che sono limitati ai miRNA noti, il sequenziamento dei miRNA (miRNA-Seq) consente una profilazione completa, ad alto rendimento e imparziale di miRNA noti e nuovi in diversi tipi di campioni e condizioni⁴. Il miRNA-Seq offre sensibilità, accuratezza e gamma dinamica superiori, rendendolo un metodo preferito per studiare i modelli di espressione dei miRNA e scoprire meccanismi di regolazione in contesti fisiologici e patologici⁵. Tuttavia, l'analisi dei dati miRNA-Seq presenta sfide computazionali specifiche, tra cui la gestione di lunghezze di lettura brevi, la rimozione di sequenze adattatrici, la distinzione tra membri della famiglia miRNA strettamente correlati e la gestione di un'elevata ridondanza nei conteggi delle letture⁶. Queste caratteristiche richiedono un flusso di lavoro analitico accuratamente progettato e standardizzato.

Sebbene siano state sviluppate varie pipeline e strumenti software per l'analisi dei dati miRNA-Seq, molti di essi si basano su interfacce utente grafiche o flussi di lavoro fissi che limitano la flessibilità e la riproducibilità⁷. Al contrario, l'ambiente di programmazione R fornisce una piattaforma potente e personalizzabile per l'analisi bioinformatica⁸. R offre un ricco ecosistema di pacchetti per la modellazione statistica, la visualizzazione dei dati e l'integrazione con database biologici. Ciò consente agli utenti di eseguire analisi complete e riproducibili in modo trasparente e basato su script. Inoltre, la natura modulare dei flussi di lavoro R consente ai ricercatori di personalizzare ogni fase in base a specifici requisiti sperimentali, dalla pre-elaborazione dei dati grezzi all'interpretazione funzionale.

In questo protocollo, presentiamo un flusso di lavoro di analisi miRNA-Seq verificato e completo implementato interamente in R, con l'obiettivo di fornire una soluzione riproducibile e adattabile all'utente per i ricercatori che lavorano con i dati di espressione dei miRNA. Il flusso di lavoro inizia con il controllo di qualità e il taglio dell'adattatore delle letture di sequenziamento grezze, seguito dall'allineamento a un genoma di riferimento o a sequenze di miRNA note. Le fasi successive includono la quantificazione dei conteggi delle letture, la normalizzazione, l'analisi dell'espressione differenziale, la previsione del gene bersaglio, l'arricchimento funzionale e la visualizzazione della rete. Il flusso di lavoro incorpora diversi pacchetti R ampiamente utilizzati e ben mantenuti, garantendo affidabilità e compatibilità con gli aggiornamenti e le estensioni futuri.

Uno dei punti di forza di questo protocollo risiede nella sua capacità di andare oltre i risultati di espressione differenziale e fornire un'interpretazione biologica significativa. Integrando database curati di interazioni miRNA-mRNA convalidate e previste, il flusso di lavoro consente agli utenti di identificare i geni bersaglio biologicamente rilevanti. Questi bersagli possono quindi essere sottoposti ad analisi di ontologia genica e di arricchimento del percorso per scoprire i processi biologici e i percorsi molecolari interessati. Nella fase finale, le reti di interazione miRNA-mRNA possono essere visualizzate utilizzando strumenti esterni come Cytoscape⁹, fornendo informazioni sul panorama normativo e identificando i miRNA hub chiave con potenziale importanza funzionale.

Questo metodo è stato applicato con successo in contesti di ricerca clinica, compresi gli studi sulle malattie renali, dove i miRNA circolanti fungono da promettenti biomarcatori per la diagnosi e la prognosi¹⁰. Tuttavia, il design modulare e flessibile del flusso di lavoro lo rende adatto a un'ampia gamma di applicazioni, tra cui la modellazione di malattie, gli studi sulla risposta ai farmaci, la biologia dello sviluppo e la genomica comparativa. I ricercatori possono facilmente adattare il flusso di lavoro per adattarsi ad annotazioni specie-specifiche, condizioni sperimentali o livelli aggiuntivi di dati omici.

Offrendo una soluzione open source basata su script, questa pipeline centrata sull'R affronta molte delle principali limitazioni associate agli strumenti miRNA-Seq esistenti, tra cui la personalizzazione limitata, la dipendenza da interfacce grafiche non trasparenti, la mancanza di supporto per organismi non modello, la scarsa riproducibilità dovuta all'assenza di controllo di versione e la difficoltà di integrazione con i framework di analisi statistica e funzionale a valle. Consente il pieno controllo dei parametri di elaborazione dei dati, incoraggia la riproducibilità attraverso il codice a controllo di versione e promuove la trasparenza nella ricerca bioinformatica. Poiché l'importanza dei miRNA continua a crescere nel contesto della biologia dei sistemi e della medicina traslazionale, avere accesso a un framework di analisi affidabile e adattabile diventa sempre più essenziale.

Protocol

NOTA: I materiali con collegamenti software sono elencati nella Tabella dei materiali.

1. Preparare campioni di RNA e librerie di sequenze

NOTA: Eseguire l'estrazione e il sequenziamento dell'RNA al di fuori di questo flusso di lavoro computazionale. Esiste più di un modo per analizzare i dati di sequenziamento dei miRNA. Questa sezione fornisce il contesto di una pratica.

1. Estrarre l'RNA totale: estrarre l'RNA totale dai campioni biologici utilizzando un kit ottimizzato per l'isolamento di piccoli RNA (ad esempio, kit di isolamento dei miRNA). Seguire attentamente il protocollo del produttore. Assicurati di utilizzare materiali di consumo privi di RNasi e di conservare i campioni sul ghiaccio per ridurre al minimo la degradazione.
2. Valutare l'integrità e la quantità dell'RNA: eseguire 1-2 μL dell'RNA estratto su un bioanalizzatore o un dispositivo equivalente. Controllare il numero di integrità dell'RNA (RIN) e assicurarsi che sia ≥ 7.0 per un sequenziamento affidabile. Registrare la concentrazione utilizzando uno spettrofotometro o un fluorimetro.
3. Costruire piccole librerie di RNA: preparare librerie di sequenziamento da 1 μg di RNA totale utilizzando un kit di preparazione di librerie di piccole dimensioni di RNA-seq commerciale. Seguire il protocollo del kit per legare gli adattatori, trascrivere inversamente e amplificare il cDNA. Purificare i prodotti della PCR con la selezione delle dimensioni (ad es. inserto da 18-30 nt) per arricchirli per i frammenti di miRNA.
4. Librerie di sequenze: caricare le librerie su una piattaforma di sequenziamento ad alta produttività. Impostare la configurazione di esecuzione per il sequenziamento single-end con una lunghezza di lettura di ~50 bp. Assicurarsi che ogni campione generi circa 10 milioni di letture non elaborate per ottenere una profondità sufficiente.
5. Esportazione dell'output del sequenziamento: Dopo il sequenziamento, esportare i dati grezzi come file FASTQ utilizzando il software di output dati dello strumento. Verificare che la directory di output contenga sia le letture della sequenza che i file del punteggio di qualità associati. Archiviare i file FASTQ in una directory strutturata per l'analisi a valle.

2. Pre-elaborare le letture grezze ed eseguire il controllo di qualità

1. Sequenze dell'adattatore di ritaglio
   1. Installa e configura Cutadapt o fastp.
   2. Eseguire il taglio dell'adattatore su ogni file FASTQ utilizzando il comando seguente:
      cutadapt -a XXXX -o trimmed_reads.fastq raw_reads.fastq
      NOTA: La sequenza di adattatori seguita dopo '-a' deve essere applicata con il kit di preparazione della libreria di sRNA specifico. '-o' definisce il nome del file di output con il file di input seguito.
2. Valutare la qualità di lettura
   1. Utilizza FastQC per generare report di controllo qualità:
      fastqc trimmed_reads.fastq
   2. Esamina i punteggi di qualità per base, la distribuzione della lunghezza di lettura e la contaminazione dell'adattatore: apri il report HTML FastQC generato per ogni file FASTQ in un browser Web. Esamina passo dopo passo i seguenti moduli:
      1. Qualità per sequenza di base: assicurarsi che la maggior parte delle basi rientri nell'area verde (punteggio Phred ≥30). Cerca qualsiasi calo di qualità all'estremità 3' che potrebbe indicare errori di sequenziamento.
      2. Leggere la distribuzione della lunghezza: confermare che la distribuzione corrisponda alla dimensione prevista dell'inserto (ad esempio, 18-30 nt per i miRNA). Verificare che non siano presenti picchi imprevisti.
      3. Contenuto dell'adattatore: verificare che le sequenze dell'adattatore siano state tagliate in modo efficace. Verificare che la percentuale di contaminazione dell'adattatore sia vicina allo zero dopo la rifilatura.
   3. Salva il report di riepilogo FastQC e contrassegna eventuali campioni con metriche di scarsa qualità per il ritaglio o l'esclusione da ulteriori analisi.

3. Mappare le letture e generare matrici di conteggio

1. Allinea le letture al riferimento
   1. Scarica il file FASTA del genoma di riferimento o delle sequenze di miRNA mature (ad esempio, da miRBase)¹¹. Esempio:
      wget ftp://mirbase.org/pub/mirbase/CURRENT/mature.fa
   2. Indicizza il genoma di riferimento usando Bowtie.
      1. Apri un terminale ed esegui il seguente comando per creare l'indice:
         bowtie-build reference.fa reference_index
      2. Sostituisci reference.fa con il nome effettivo del file FASTA.
      3. Sostituire reference_index con il prefisso desiderato per l'indice.
      4. Assicurati che Bowtie generi più file di indice (ad esempio, .ebwt). Verificare che questi file siano presenti nella directory di lavoro, in quanto sono necessari per l'allineamento.
   3. Allinea le letture utilizzando Bowtie con parametri adatti per letture brevi. Esempio:
      bowtie -v 0 -a --best --strata reference_index trimmed_reads.fastq > aligned_reads.sam
      NOTA: Il file di input è trimmed_reads.fastq e il file di output è aligned_reads.sam. '-v 0' significa che non permettiamo errori nell'intera lettura. '-a -best -strata' significa che buttiamo via qualsiasi allineamento che abbia più discrepanze del migliore.
2. Quantificare l'espressione dei miRNA
   1. Converti i file SAM in formato BAM utilizzando SAMtools.
      samtools view -S -b aligned_reads.sam > aligned_reads.bam
      NOTA: Il file di input aligned_reads.sam è il risultato dell'ultimo comando. Il file di output aligned_reads.bam viene preparato per l'analisi successiva.
      1. Utilizzare SAMtools per comprimere e ordinare il file di allineamento:
         samtools sort aligned_reads.bam -o aligned_reads_sorted.bam
samtools index aligned_reads_sorted.bam
      2. Assicurarsi che il primo comando converta il file SAM in formato BAM.
      3. Assicurarsi che il secondo comando ordini il file BAM in base alle coordinate genomiche.
      4. Assicurarsi che il terzo comando generi un file di indice (.bai), necessario per le analisi a valle.
      5. Verificare che sia il file BAM ordinato che il relativo indice siano stati creati correttamente prima di procedere alla quantificazione.
   2. Utilizzare featureCounts o HTSeq-count per generare una matrice di conteggio utilizzando l'annotazione miRNA GTF:
      featureCounts -a miRNA.gtf -o counts.txt aligned_reads.bam
      NOTA: featureCounts quantifica le letture nel file di input aligned_reads.bam in base a miRNA.gtf e restituisce counts.txt.

4. Eseguire l'analisi dell'espressione differenziale in R

1. Dati di conteggio del carico
   1. Importare la matrice di conteggio e i metadati del campione in 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)
      NOTA: I dati di conteggio (counts.csv) e il gruppo di campioni (metadata.csv) devono essere forniti a DESeq2. La "condizione" qui chiarisce il gruppo di campioni forniti. Per le esigenze individuali, fare riferimento al manuale di DESeq2¹².
2. Normalizzare e trasformare i dati
   1. Normalizzare i dati di conteggio utilizzando il metodo predefinito di DESeq2:
      dds <- DESeq(dds)
   2. Eseguire la trasformazione stabilizzante della varianza:
      vsd <- vst(dds, blind=FALSE)
   3. Visualizzare il clustering dei campioni utilizzando PCA:
      plotPCA(vsd, intgroup="condition")
3. Identificare i miRNA differenzialmente espressi
   1. Estrarre e ordinare i risultati delle espressioni differenziali:
      res <- results(dds)
resOrdered <- res[order(res$pvalue), ]
      NOTA: Riordiniamo il file dei risultati 'res' in base al valore di pvalue.
      summary(res)
   2. Filtro per miRNA espressi in modo significativamente differenziale (valore p < 0,05, |log2FC| > 1).
      sig_miRNA <- subset(res, pvalue < 0.05 & abs(log2FC) > 1)
      NOTA: Esistono più soglie per filtrare i miRNA significativamente modificati. I requisiti 'p-value < 0.05, |log2FC| > 1' sono ampiamente applicati. Le soglie possono essere regolate per i singoli dati.
4. Visualizzare le modifiche alle espressioni
   1. Installa e carica il pacchetto EnhancedVolcano.
   2. Crea la trama del vulcano:
      library(EnhancedVolcano)
EnhancedVolcano(res,
lab = rownames(res),
x = 'log2FoldChange',
y = 'pvalue',
title = 'Differentially Expressed miRNAs')

5. Predire i geni bersaglio dei miRNA

1. Interrogare i database
   1. Utilizza risorse online come TargetScan, miRDB e miRTarBase per cercare microRNA specifici e recuperare geni bersaglio.
   2. Concentrati su obiettivi convalidati sperimentalmente per una maggiore sicurezza.
2. Automatizza la previsione in R
   1. Caricare il pacchetto multiMiR ed eseguire query sulle destinazioni convalidate:
      library(multiMiR)
target_results <- get_multimir(mirna = c("hsa-miR-21-5p"), table = "validated")
      NOTA: Qui prendiamo ad esempio "has-miR-21-5p" e otteniamo bersagli convalidati.
   2. Estrarre simboli unici del gene bersaglio per l'analisi dell'arricchimento:
      genes <- unique(target_results@data$target_symbol)

6. Condurre un'analisi dell'arricchimento funzionale

1. Eseguire l'arricchimento GO
   1. Strumenti di arricchimento del carico:
      library(clusterProfiler)
library(org.Hs.eg.db)
      NOTA: Qui carichiamo un database contenente annotazioni sul genoma umano, utili per convertire gli identificatori genici comuni.
   2. Esegui l'analisi dell'arricchimento GO per i processi biologici:
      ego <- enrichGO(gene = genes,
OrgDb = org.Hs.eg.db,
keyType = "SYMBOL",
ont = "BP",
pAdjustMethod = "BH",
pvalueCutoff = 0.05)
dotplot(ego)
      NOTA: Per l'uso di enrichGO, forniamo un elenco di geni e chiariamo che sono nominati come "SIMBOLO" qui. Abbiamo condotto l'arricchimento del processo biologico, che corrisponde al parametro 'ont = "BP" '. Eseguiamo più correzioni di test, quindi designiamo pAdjustMethod = "BH". Per una soglia significativa, scegliamo pvalueCutoff = 0.05.dotplot presenta un risultato visualizzato. Per altre opzioni individuali, fare riferimento al manuale di clusterProfiler¹³.
2. Eseguire l'arricchimento del percorso KEGG
   1. Eseguire l'arricchimento KEGG:
      ekegg <- enrichKEGG(gene = genes, organism = 'hsa')
dotplot(ekegg)
      NOTA: Per l'utilizzo di enrichKEGG, forniamo un elenco di geni e chiariamo l'organismo come umano ("ha"). dotplot presenta un risultato visualizzato. Per altre opzioni individuali, fare riferimento al manuale di clusterProfiler¹³.

7. Costruire e visualizzare la rete di interazione miRNA - mRNA

1. Esportare i dati per la visualizzazione della rete
2. Crea un frame di dati di coppie di geni miRNA-bersaglio in base ai bersagli generati da TargetScan, miRDB o miRTarBase.
3. Scrivi la tabella di rete in CSV:
   write.csv(miRNA_target_pairs, "network.csv")
4. Importazione in Cytoscape
   1. Apri Cytoscape e importa la tabella di rete.
   2. Visualizza la rete utilizzando il layout circolare o a forza diretta.
   3. Analizzare le proprietà topologiche (ad esempio, il grado di centralità) per identificare i miRNA hub.

Results

Abbiamo scaricato la matrice di espressione dei microRNA da GSE133530 e condotto direttamente l'analisi dell'espressione differenziale. È stato fornito un esempio di script R analitico per il set di dati nel file supplementare 1. Il set di dati ha eseguito la profilazione globale dei miRNA su 16 cisti renali di diverse dimensioni (cisti minime: meno di 1-5 mL, n = 10; cisti medie: tra 10-25 mL, n = 4; cisti grandi: maggiori di 50 mL, n = 4) e tessuto minimamente cistico (...

Discussion

L'analisi dei dati miRNA-Seq presenta sfide distinte a causa delle dimensioni ridotte e della ridondanza delle letture, rendendo fondamentale un rigoroso controllo di qualità e la pre-elaborazione. Uno dei passaggi più importanti del flusso di lavoro è la rifilatura dell'adattatore. Poiché i miRNA sono lunghi circa 22 nucleotidi, le sequenze di adattatori possono facilmente dominare le letture se non vengono rimosse correttamente. La mancata esecuzione di un taglio accurato può causare d...

Materials

List of materials used in this article

Name	Company	Catalog Number	Comments
Agilent-021827 Microarray di miRNA umano	Agilent	/	Un array commerciale per la profilazione di microRNA di campioni umani
Cravatta a farfalla	Università Johns Hopkins	http://bowtie-bio.sourceforge.net/index.shtml	Uno strumento software per l'allineamento delle letture di sequenziamento a lunghe sequenze di riferimento
clusterProfiler (pacchetto R)	Bioconduttore	https://bioconductor.org/packages/clusterProfiler/	Un pacchetto R progettato per l'analisi dell'arricchimento funzionale e la visualizzazione di dati biologici ad alto rendimento.
Cutadapt	Open Source	https://cutadapt.readthedocs.io	Uno strumento da riga di comando che rimuove sequenze di adattatori, primer, code poly-A e altri frammenti indesiderati dalle letture di sequenziamento ad alto rendimento.
Cytoscape	Consorzio Cytoscape	https://cytoscape.org/	Una piattaforma software open-source progettata per la visualizzazione e l'analisi di reti biologiche complesse.
DESeq2 (pacchetto R)	Bioconduttore	https://bioconductor.org/packages/DESeq2/	Un pacchetto R progettato per l'analisi dell'espressione genica differenziale dei dati di conteggio
EnhancedVolcano (pacchetto R)	Bioconduttore	https://bioconductor.org/packages/EnhancedVolcano/	Un pacchetto R progettato per creare grafici vulcanici di qualità da pubblicazione.
Controllo rapido	Babraham Bioinformatica	https://www.bioinformatics.babraham.ac.uk/projects/fastqc/	Uno strumento di controllo qualità open source per dati di sequenziamento ad alto rendimento.
featureCounts	Lettura secondaria / SourceForge	http://subread.sourceforge.net/	Un programma utilizzato per il conteggio delle letture mappate alle caratteristiche genomiche
HTSeq-count	Pacchetto Python	https://htseq.readthedocs.io	Uno strumento da riga di comando che conta il numero di letture di sequenziamento ad alto rendimento allineate che si sovrappongono a caratteristiche genomiche come geni o esoni. Io
Illumina Human v2 MicroRNA expression beadchip	Illumina	/	Un array commerciale per la profilazione di microRNA di campioni umani
multiMiR (pacchetto R)	Bioconduttore	https://bioconductor.org/packages/multiMiR/	Un pacchetto R che fornisce la più ampia raccolta integrata di microRNA previsti e convalidati sperimentalmente. Interazioni mirate insieme alle loro associazioni a malattie e farmaci.
Org. Hs.eg.db (pacchetto R)	Bioconduttore	https://bioconductor.org/packages/org.Hs.eg.db/	Un pacchetto di annotazioni progettato per la ricerca genomica umana (Homo sapiens).
R Software	Progetto R	https://www.r-project.org/	Un progetto open-source per il calcolo statistico
Rstudio	Poni PBC	/	Un ambiente di sviluppo integrato aiuta a essere più produttivi con R e Python
Strumenti SAM	Open Source	http://www.htslib.org/	Un pacchetto software per la manipolazione dei dati di sequenziamento di nuova generazione (NGS).