In Silico Identifizierung und Charakterisierung von circRNAs während Wirt-Pathogen-Interaktionen

Das hier eingereichte Protokoll erläutert die vollständige In-silico-Pipeline , die für die Vorhersage und funktionelle Charakterisierung von circRNAs aus RNA-Sequenzierungs-Transkriptomdaten zur Untersuchung von Wirt-Pathogen-Interaktionen erforderlich ist.

Zirkuläre RNAs (circRNAs) sind eine Klasse nicht-kodierender RNAs, die durch Back-Sppleißen gebildet werden. Diese circRNAs werden vor allem auf ihre Rolle als Regulatoren verschiedener biologischer Prozesse untersucht. Bemerkenswert ist, dass neue Erkenntnisse zeigen, dass Wirts-circRNAs bei einer Infektion mit Krankheitserregern (z. B. Influenza und Coronaviren) differentiell exprimiert werden können, was auf eine Rolle von circRNAs bei der Regulierung der angeborenen Immunantwort des Wirts hindeutet. Untersuchungen zur Rolle von circRNAs bei pathogenen Infektionen sind jedoch begrenzt durch die Kenntnisse und Fähigkeiten, die erforderlich sind, um die notwendigen bioinformatischen Analysen durchzuführen, um DE-circRNAs aus RNA-Sequenzierungsdaten (RNA-seq) zu identifizieren. Die bioinformatische Vorhersage und Identifizierung von circRNAs ist von entscheidender Bedeutung vor jeder Verifizierung und funktionellen Studien mit kostspieligen und zeitaufwändigen Nasslabortechniken. Um dieses Problem zu lösen, wird in diesem Manuskript ein Schritt-für-Schritt-Protokoll zur In-silico-Vorhersage und Charakterisierung von circRNAs unter Verwendung von RNA-seq-Daten bereitgestellt. Das Protokoll kann in vier Schritte unterteilt werden: 1) Vorhersage und Quantifizierung von DE-circRNAs über die CIRIquant-Pipeline; 2) Annotation mittels circBase und Charakterisierung von DE circRNAs; 3) Vorhersage der CircRNA-miRNA-Interaktion durch die Circr-Pipeline; 4) Analyse der funktionellen Anreicherung von circRNA-Elterngenen unter Verwendung von Gene Ontology (GO) und Kyoto Encyclopedia of Genes and Genomes (KEGG). Diese Pipeline wird nützlich sein, um die zukünftige In-vitro - und In-vivo-Forschung voranzutreiben, um die Rolle von circRNAs bei Wirt-Pathogen-Interaktionen weiter zu entschlüsseln.

Wirt-Pathogen-Interaktionen stellen ein komplexes Zusammenspiel zwischen den Krankheitserregern und den Wirtsorganismen dar, das die angeborenen Immunantworten des Wirts auslöst, die schließlich zur Entfernung eindringender Krankheitserreger führen ^1,2. Bei pathogenen Infektionen wird eine Vielzahl der Immungene des Wirts reguliert, um die Vermehrung und Freisetzung von Krankheitserregern zu hemmen. Zu den häufigen Interferon-stimulierten Genen (ISGs), die bei pathogenen Infektionen reguliert werden, gehören ADAR1, IFIT1, IFIT2, IFIT3, ISG20, RIG-I und OASL ^3,4. Studien haben gezeigt, dass neben proteinkodierenden Genen auch nicht-kodierende RNAs wie lange nicht-kodierende RNAs (lncRNAs), microRNAs (miRNAs) und zirkuläre RNAs (circRNAs) eine Rolle spielen und gleichzeitig bei pathogenen Infektionen reguliert werden ^5,6,7. Im Gegensatz zu proteinkodierenden Genen, die Proteine hauptsächlich als funktionelle Moleküle kodieren, sind nicht-kodierende RNAs (ncRNAs) als Regulatoren von Genen auf transkriptioneller und posttranskriptioneller Ebene bekannt. Allerdings sind Studien, in denen nicht-kodierende RNAs, insbesondere circRNAs, an der Regulation der Immungene des Wirts beteiligt sind, im Vergleich zu den proteinkodierenden Genen nicht gut berichtet.

CircRNAs zeichnen sich weitgehend durch ihre kovalent geschlossene Endlosschleifenstruktur aus, die durch einen nicht-kanonischen Spleißprozess namens Back-Spleißen⁸ erzeugt wird. Im Gegensatz zum Spleißen verwandter linearer RNAs wird beim Back-Splicing die nachgeschaltete Donorstelle an die stromaufwärts gelegene Akzeptorstelle ligiert, wodurch eine kreisförmige Struktur entsteht. Derzeit werden drei verschiedene Back-Spleiß-Mechanismen für die Biogenese von circRNAs vorgeschlagen. Dabei handelt es sich um RNA-bindende Proteine (RBP)-vermittelte Zirkularisierung 9,10, Intron-Paarungs-getriebene Zirkularisierung 11 und Lariat-getriebene Zirkularisierung^12,13,14. Da circRNAs in einer kreisförmigen Struktur Ende-zu-Ende verbunden sind, neigen sie dazu, von Natur aus resistent gegen normale Exonuklease-Verdauungen zu sein und gelten daher als stabiler als ihre linearen Gegenstücke¹⁵. Ein weiteres gemeinsames Merkmal von circRNAs ist die zell- oder gewebetypspezifische Expression in Wirten¹⁶.

Wie ihre einzigartige Struktur und ihre zell- oder gewebespezifische Expression vermuten lassen, haben circRNAs wichtige biologische Funktionen in Zellen übernommen. Bis heute ist eine der herausragenden Funktionen von circRNAs ihre Rolle als microRNA (miRNA)-Schwämme^17,18. Diese regulatorische Rolle von circRNAs erfolgt durch die komplementäre Bindung von circRNA-Nukleotiden an die Seed-Region von miRNAs. Eine solche circRNA-miRNA-Interaktion hemmt die normalen regulatorischen Funktionen der miRNAs auf den Ziel-mRNAs und reguliert so die Expression der Gene^19,20. Darüber hinaus ist bekannt, dass circRNAs die Genexpression regulieren, indem sie mit RNA-bindenden Proteinen (RBPs) interagieren und RNA-Protein-Komplexe bilden²¹. Obwohl circRNAs als nicht-kodierende RNAs klassifiziert werden, gibt es auch Hinweise darauf, dass circRNAs als Vorlagen für die Proteintranslation fungieren können^22,23,24.

In jüngster Zeit wurde gezeigt, dass circRNAs eine zentrale Rolle bei der Regulierung der Wirt-Pathogen-Interaktionen spielen, insbesondere zwischen Wirt und Viren. Im Allgemeinen wird angenommen, dass Wirts-circRNAs bei der Regulierung der Immunantwort des Wirts helfen, um die eindringenden Krankheitserreger zu eliminieren. Ein Beispiel für circRNA, die die Immunantwort des Wirts fördert, ist circRNA_0082633, berichtet von Guo et ^al.25. Diese circRNA verstärkt die Signalübertragung von Typ-I-Interferon (IFN) in A549-Zellen, was dazu beiträgt, die Replikation von Influenzaviren zu unterdrücken²⁵. Darüber hinaus berichteten Qu et al. auch über eine humane intronische circRNA, genannt circRNA AIVR, die die Immunität fördert, indem sie die Expression des CREB-bindenden Proteins (CREBBP), einem Signalwandler von IFN-β^{reguliert 26,27}. Es gibt jedoch auch circRNAs, von denen bekannt ist, dass sie die Pathogenese von Krankheiten bei einer Infektion fördern. So berichteten Yu et al. kürzlich über die Rolle einer circRNA, die aus der GATA-Zinkfingerdomäne gespleißt wird, die das 2A-Gen (circGATAD2A) enthält, bei der Förderung der Replikation des H1N1-Virus durch die Hemmung der Autophagie der Wirtszelle²⁸.

Um circRNAs effektiv untersuchen zu können, wird in der Regel ein genomweiter circRNA-Vorhersagealgorithmus implementiert, gefolgt von einer In-silico-Charakterisierung der vorhergesagten circRNA-Kandidaten, bevor funktionelle Studien durchgeführt werden können. Ein solcher bioinformatischer Ansatz zur Vorhersage und Charakterisierung von circRNAs ist kostengünstiger und zeiteffizienter. Es hilft, die Anzahl der Kandidaten zu verfeinern, die funktionell untersucht werden sollen, und könnte möglicherweise zu neuen Erkenntnissen führen. In dieser Arbeit stellen wir ein detailliertes bioinformatisches Protokoll für die in silico Identifizierung, Charakterisierung und funktionelle Annotation von circRNAs während der Wirt-Pathogen-Interaktionen zur Verfügung. Das Protokoll umfasst die Identifizierung und Quantifizierung von circRNAs aus RNA-Sequenzierungsdatensätzen, die Annotation über circBase und die Charakterisierung der circRNA-Kandidaten in Bezug auf circRNA-Typen, Anzahl überlappender Gene und vorhergesagte circRNA-miRNA-Interaktionen. Diese Studie liefert auch die funktionelle Annotation der circRNA-Elterngene durch Genontologie (GO) und die Anreicherungsanalyse der Kyoto Encyclopedia of Genes and Genomes (KEGG).

In diesem Protokoll wurden de-identifizierte ribosomale RNA (rRNA)-depletierte RNA-seq-Bibliotheksdatensätze, die aus den mit dem Influenza-A-Virus infizierten menschlichen Makrophagenzellen erstellt wurden, aus der Gene Expression Omnibus (GEO)-Datenbank heruntergeladen und verwendet. Die gesamte bioinformatische Pipeline von der Vorhersage bis zur funktionellen Charakterisierung von circRNAs ist in Abbildung 1 zusammengefasst. Jeder Teil der Pipeline wird in den folgenden Abschnitten näher erläutert.

1. Vorbereitung, Download und Einrichtung vor der Datenanalyse

HINWEIS: Alle in dieser Studie verwendeten Softwarepakete sind kostenlos und quelloffen.

1. Herunterladen der erforderlichen Tools auf der Linux-Plattform
   1. Laden Sie die erforderliche Software und die Tools, die im Materialverzeichnis aufgeführt sind, herunter und installieren Sie sie auf einem Linux-Hochleistungscomputer, indem Sie die Anweisungen des Entwicklers befolgen.
      HINWEIS: Die meisten Tools und Software verfügen über eigene Online-GitHub-Seiten oder Dokumentationen mit Anweisungen zur Installation und Verwendung ihrer Tools (siehe Materialtabelle).
   2. Laden Sie die gewünschten RNA-seq-Datensätze für die circRNA-Detektion und -Analyse von Sequenzarchiv-Websites herunter (z. B. European Nucleotides Archive und Gene Expression Omnibus).
   3. Laden Sie die Referenzgenome (FASTA-Format) und Annotationsdateien (GTF/GFF3-Format) herunter, die mit dem Host kompatibel sind, von dem der RNA-seq-Datensatz erstellt wurde. Referenzgenome und Annotationsdateien des Wirts sind in der Regel in Online-Genombrowsern wie dem National Center for Biotechnology Information (NCBI), der University of California Santa Cruz (UCSC) und den Ensembl-Websites zu finden.
2. Qualitätsprüfung von RNA-seq
   1. Geben Sie die FASTQ-Dateien in das FASTQC-Programm ein, um die Qualität der RNA-Sequenzen zu bestimmen. Wenn die Qualität der FASTQ-Dateien niedrig ist (z. B. 29,30 erforderlich sein.

2. Vorhersage und differentielle Expressionsanalyse von circRNAs mit CIRIquant

HINWEIS: Ein ausführlicheres Handbuch zur Installation und Durchführung der Differentialexpressionsanalyse finden Sie im Abschnitt Codeverfügbarkeit des CIRIquant-Papiers³¹. Die ergänzenden Daten enthalten auch einige der grundlegenden Befehle, die in diesem Protokoll verwendet werden.

1. CircRNA-Vorhersagen
   1. Indizieren Sie zuerst das Referenzgenom des Wirts mit BWA- und HISAT2-Alignern. Führen Sie dann auf einem Linux-Terminal die Befehle bwa index 32 und hisat2-build³³ im Verzeichnis des Referenzgenoms des Hosts aus, um es zu indizieren.
   2. Bereiten Sie als Nächstes eine YML-Konfigurationsdatei vor, die den Namen der Datei, den Pfad der Tools (BWA, HISAT2, stringtie³⁴, samtools³⁵), den Pfad zu den heruntergeladenen Referenzdateien (FASTA-Dateien des Referenzgenoms des Hosts, Annotationsdateien) und den Pfad zu den Indexdateien aus Schritt 2.1.1 enthält.
   3. Führen Sie das CIRIquant-Tool vom Terminal aus mit den Standard- oder manuellen Parametern aus. Der Benutzer kann den Bibliothekstyp (entweder stranged oder non-stranged) der RNA-seq-Daten angeben, wenn er das CIRIquant-Tool ausführt.
      HINWEIS: Der Bibliothekstyp der RNA-seq-Daten kann bestimmt werden, indem der Typ des verwendeten Bibliotheksvorbereitungskits bekannt ist. Wenn die Identität des Bibliothekspräparationskits unbekannt ist, kann ein bioinformatisches RNA-seq-Kontrollpaket namens RSeQC³⁶ verwendet werden, um die Strandedness von RNA-seq-Daten zu bestimmen.
2. Differentielle Expressionsanalyse
   HINWEIS: Das CIRIquant-Paket enthält prep_CIRIquant, prepDE.py und CIRI_DE_replicate; Daher sind für diese drei Tools keine zusätzlichen Downloads erforderlich.
   1. Bereiten Sie eine Textdatei (.lst) mit einer Liste von Daten vor, die Folgendes enthält:
      1^{. Spalte} : IDs der in Schritt 2.1.3 verwendeten RNA-seq-Daten
      2^. Spalte: Pfad zu den GTF-Dateien, die von CIRIquant ausgegeben werden
      3. Spalte: Gruppierung der RNA-seq-Daten, unabhängig davon, ob es sich um eine Kontroll- oder eine behandelte Gruppe handelt^.
   2. Ein Beispiel finden Sie in Tabelle 1 unten.
      HINWEIS: Es ist nicht notwendig, die Kopfzeilen einzufügen, da sie nur als Referenz dienen.
   3. Führen Sie auf dem Linux-Terminal prep_CIRIquant mit der in Schritt 2.2.1 vorbereiteten Textdatei (.lst) als Eingabe aus. Bei der Ausführung wird eine Liste von Dateien generiert: library_info.csv, circRNA_info.csv, circRNA_bsj.csv und circRNA_ratio.csv.
   4. Bereiten Sie eine zweite Textdatei mit einer Liste von Daten vor, die die RNA-seq-IDs und den Pfad zu ihrer jeweiligen StringTie-Ausgabe enthält. Das Dateilayout muss der Textdatei in Schritt 2.2.1 ohne die Gruppierungsspalte ähneln.
   5. Führen Sie prepDE.py mit der in Schritt 2.2.4 vorbereiteten Textdatei als Eingabe aus, um die Genzählmatrixdateien zu generieren.
   6. Führen Sie CIRI_DE_replicate mit den library_info.csv - und circRNA_bsj.csv Dateien aus Schritt 2.2.3 und der gene_count_matrix.csv Datei aus Schritt 2.2.5 als Eingaben aus, um die endgültige circRNA_de.tsv-Datei auszugeben.
3. Filterung von DE-circRNAs
   1. Verwenden Sie R (im Computerterminal oder RStudio) oder eine beliebige Tabellenkalkulationssoftware (z. B. Microsoft Excel), um die aus Schritt 2.2.6 generierte Datei circRNA_de.tsv zu öffnen, um die Anzahl der differentiell exprimierten (DE) circRNAs zu filtern und zu bestimmen.
   2. Filtern Sie die DE-circRNAs nach den Kriterien LogFC > |2| und FDR < 0,05.
   3. Erstellen Sie eine Datei mit dem Namen DE_circRNAs.txt , um die Informationen von DE-circRNAs zu speichern.

3. Charakterisierung und Annotation von vorhergesagten DE-circRNAs

1. Annotationsstatus von DE-circRNAs
   1. Laden Sie die Datei mit dem Namen DE_circRNAs.txt in RStudio , die aus der Liste der DE-circRNAs besteht, die aus Schritt 2.3.3 gefiltert wurden. Fügen Sie weitere Informationen hinzu, wie z. B. die genomischen Positionen (Chr, Start, Ende), die Strangorientierung (+ oder -), den Gennamen und den circRNA-Typ. Bevor Sie fortfahren, konvertieren Sie die genomischen Startkoordinaten der circRNA von CIRIquant in 0-basiert, indem Sie 1 Basenpaar subtrahieren.
      HINWEIS: Die anderen oben genannten Informationen können den GTF-Dateien entnommen werden, die von CIRIquant (Supplementary File 1) ausgegeben werden.
   2. Bestimmen Sie den Annotationsstatus der vorhergesagten DE-circRNAs, indem Sie eine Bibliothek herunterladen, die die genomischen Positionen der in der circRNA-Datenbank (z.B. circBase) hinterlegten circRNAs enthält.
      HINWEIS: Stellen Sie sicher, dass die Genomversion, die zur Vorhersage der circRNAs verwendet wird, mit der circRNA-Datenbankbibliothek identisch ist, bevor Sie den Vergleich durchführen. Die hier verwendete circBase-Datendatei ist in dem in Github (https://github.com/bicciatolab/Circr)³⁷ bereitgestellten Laufwerksordner frei verfügbar.
   3. Sobald die beiden Dateien aus Schritt 3.1.1 und Schritt 3.1.2 vorbereitet sind, führen Sie das in Ergänzende Datei 1 angegebene R-Skript aus. Chromosomale Positionen von DE-circRNAs werden an die Bibliothek abgefragt, bevor ihnen der Status Annotated oder Unannotated zugewiesen wird.
2. Charakterisierung von DE-circRNAs
   1. Verwenden Sie R und andere Tabellenkalkulationsprogramme, um die Anzahl der circRNAs entsprechend den circRNA-Typen (d. h. Exon, Intron, intergen und Antisense) und die Anzahl der Gene, über die sich die circRNAs erstrecken (1 oder >1), zusammenzufassen (Ergänzungsdatei 1).HINWEIS: CIRIquant kann nur vier Arten von circRNAs nachweisen (Exon, Intron, intergen und Antisense). Exon-Intron-circRNAs, auch ElciRNAs genannt, können von CIRIquant nicht nachgewiesen werden.

4. Vorhersage der circRNA-miRNA-Interaktion mit Circr

HINWEIS: Eine ausführlichere Anleitung zur Installation und Verwendung von Circr für die circRNA-miRNA-Interaktionsanalyse finden Sie unter: https://github.com/bicciatolab/Circr³⁷.

1. Vorbereitung von Dateien
   1. Entpacken und extrahieren Sie den Inhalt der Circr.zip Datei, nachdem Sie sie von der Circr GitHub-Seite mit der entsprechenden Software wie "WinRar" oder "7-zip" heruntergeladen haben, in ein neues Verzeichnis, in dem die Analyse durchgeführt wird.
   2. Installieren Sie die erforderlichen Softwareanwendungen (miRanda, RNAhybrid, Pybedtools und samtools), bevor Sie die circRNA-miRNA-Analyse durchführen.
   3. Referenzgenome und Annotationsdateien für verschiedene Organismen von Interesse, eine rRNA-Koordinatendatei, eine validierte miRNA-Interaktionsdatei und circBase-circRNA-Dateien werden vom Circr-Autor auf der Github-Seite (https://github.com/bicciatolab/Circr)³⁷ bereitgestellt. Wenn Sie auf die Unterstützungsdateien im Laufwerksordner klicken, wählen Sie den Ordner für den gewünschten Organismus, den miRNA-Ordner und die circBase-Textdatei aus und laden Sie sie herunter.
   4. Nachdem Sie die erforderlichen Dateien in Schritt 4.1.3 heruntergeladen haben, erstellen Sie ein neues Verzeichnis mit dem Namen support_files in dem in Schritt 4.1.1 genannten Verzeichnis. Entpacken Sie dann den Inhalt und extrahieren Sie ihn in das support_files Verzeichnis.
   5. Indizieren Sie die Referenzgenomdatei des interessierenden Organismus mit dem Befehl samtools faidx (Supplementary File 1).
   6. Bereiten Sie eine Eingabedatei vor, die aus den Koordinaten der gewünschten DE-circRNAs in einer tabulatorgetrennten BED-Datei besteht, wie in Tabelle 2 dargestellt.
      HINWEIS: Da circRNAs, die von CIRIquant vorhergesagt werden, nicht 0-basiert sind, ist es notwendig, 1 bp an der Startkoordinate zu minuieren (wie in Schritt 3.1.1 erwähnt), bevor sie in das BED-Format konvertiert werden. Die in Tabelle 2 gezeigten Header dienen nur als Referenz und werden in den BED-Dateien nicht benötigt.
   7. Stellen Sie an dieser Stelle sicher, dass die erwartete Ordnerstruktur für die Circr-Analyse mit der in Abbildung 2 dargestellten Struktur übereinstimmt.
2. Ausführen Circr.py
   1. Führen Sie Circr.py mit Python 3 aus, und geben Sie als Argumente die circRNA-Eingabedatei, das FASTA-Genom des gewünschten Organismus, die Genomversion des ausgewählten Organismus, die Anzahl der Threads und den Namen der Ausgabedatei in der Befehlszeile an.
   2. Wenn der betreffende Organismus nicht in dem in Schritt 4.1.3 aufgeführten Laufwerksordner bereitgestellt wird oder wenn der Benutzer es vorzieht, einen benutzerdefinierten Satz von Dateien zum Ausführen der Analyse zu haben, müssen beim Ausführen von Circr.py zusätzliche Befehle enthalten sein, die den Speicherort dieser Dateien angeben.
   3. Nachdem die Circr-Analyse abgeschlossen ist, gibt das Programm eine circRNA-miRNA-Interaktionsdatei im csv-Format aus.
   4. Filtern Sie die Ergebnisse der circRNA-miRNA-Interaktion entsprechend der benutzerspezifischen Präferenz. Für diese Studie werden die Vorhersagen mit Rstudio nach den folgenden Kriterien gefiltert:
      -Wird von allen drei Software-Tools erkannt
      -Zwei oder mehr Bindungsstellen, die sowohl von Targetscan als auch von miRanda gemeldet werden
      -Wird entweder in den Spalten "AGO" oder "validiert" identifiziert
      -Filtern Sie keine Seed-Region-Interaktionen heraus
   5. Schreiben Sie die circRNAs, die die gefilterten Bedingungen aus Schritt 4.2.3 übergeben, in eine neue Textdatei mit dem Namen circRNA_miRNA.txt. Eine solche Filterung kann die Zuverlässigkeit der vorhergesagten Wechselwirkungen erhöhen.

5. Aufbau des ceRNA-Netzwerks

HINWEIS: Eine ausführliche Anleitung zur Verwendung von Cytoscape finden Sie unter: http://manual.cytoscape.org/en/stable/ und https://github.com/cytoscape/cytoscape-tutorials/wiki#introduction

1. Download und Vorbereitung
   1. Laden Sie die neueste Version von Cytoscape³⁸ herunter von: https://cytoscape.org/download.html.
   2. Führen Sie den in Schritt 5.1.1 heruntergeladenen Installationsassistenten aus und wählen Sie den Speicherort für die Cytoscape-Software aus.
   3. Bereiten Sie eine tabulatorgetrennte Datei vor, die die interessierenden circRNAs und ihre Ziel-miRNA enthält. Die erste Spalte besteht aus dem circRNA-Namen; die zweite Spalte gibt den RNA-Typ aus der ersten Spalte an. die dritte Spalte ist die Ziel-miRNA; und die vierte Spalte gibt die Art der RNA aus der dritten Spalte an. Ein Beispiel für die Datei ist in Tabelle 3 dargestellt.
2. Aufbau der ceRNA-Netzwerkkarte
   1. Öffnen Sie die in Schritt 5.1.2 installierte Cytoscape-Software.
   2. Navigieren Sie in Cytoscape zu Datei > importieren Sie > Netzwerk aus Datei. Wählen Sie die Datei aus, die in Schritt 5.1.3 vorbereitet wurde.
   3. Wählen Sie auf der neuen Registerkarte die erste und zweite Spalte als "Quellknoten" und "Quellknotenattribut" aus, während Sie die dritte und vierte Spalte als "Zielknoten" bzw. "Zielknotenattribut" auswählen. Klicken Sie auf OK und das Netzwerk wird oben rechts in Cytoscape angezeigt.
   4. Um den visuellen Stil des Netzwerks zu ändern, klicken Sie auf die Schaltfläche Stil auf der linken Seite von Cytoscape.
   5. Drücken Sie auf den Pfeil auf der rechten Seite von Füllfarbe. Wählen Sie Typ für die Spalte und Diskrete Zuordnung für den Zuordnungstyp. Wählen Sie dann die gewünschte Farbe für jeden der RNA-Typen aus.
   6. Nachdem Sie die Farbe geändert haben, ändern Sie die Form der Knoten, indem Sie zu Form navigieren und Schritt 5.2.5 ausführen.

6. Analyse der funktionalen Anreicherung

1. Genontologie (GO) und Kyoto Encyclopedia of Genes and Genomes (KEGG) Analyse für das Elterngen der circRNAs
   1. Stellen Sie sicher, dass clusterProfiler^39,40 und org. Hs.B.db⁴¹ Pakete wurden in Rstudio installiert. Die org. Hs.eg.db^41-Paket ist ein genomweites Annotationspaket nur für Menschen. Wenn es sich bei dem interessierenden Organismus um eine andere Art handelt, siehe: https://bioconductor.org/packages/release/BiocViews.html#OrgDb
   2. Importieren Sie die DE_circRNA Informationen aus Schritt 2.3.1 in den Rstudio-Arbeitsbereich.
   3. Verwenden Sie das elterliche Gen der circRNAs, die in dieser Datei bereitgestellt werden, für die Anreicherungsanalyse in den nächsten Schritten. Wenn der Benutzer jedoch das Gensymbol in andere Formate, wie z. B. die Entrez-ID, konvertieren möchte, verwenden Sie eine Funktion wie "bitr".
   4. Wenn Sie die Gen-ID als Eingabe verwenden, führen Sie die GO-Anreicherungsanalyse mit der enrichGO-Funktion innerhalb des clusterProfiler 39,40-Pakets mit Standardparametern aus.
   5. Wenn Sie die Gen-ID als Eingabe verwenden, führen Sie die KEGG-Anreicherungsanalyse mit der enrichKEGG-Funktion innerhalb des clusterProfiler 39,40-Pakets mit den Standardparametern aus.

Das im vorherigen Abschnitt aufgeführte Protokoll wurde modifiziert und konfiguriert, um dem Linux-Betriebssystem gerecht zu werden. Der Hauptgrund dafür ist, dass die meisten Modulbibliotheken und Pakete, die an der Analyse von circRNAs beteiligt sind, nur auf der Linux-Plattform funktionieren können. In dieser Analyse wurden de-identifizierte ribosomale RNA (rRNA)-depletierte RNA-seq-Bibliotheksdatensätze, die aus den mit dem Influenza-A-Virus infizierten menschlichen Makrophagenzellen hergestellt wurden, aus der GEO-Datenbank⁴² heruntergeladen und zur Generierung der repräsentativen Ergebnisse verwendet.

CircRNA-Vorhersage und -Quantifizierung
In dieser Analyse wurden ribosomale RNA (rRNA)-depletierte RNA-seq-Bibliotheksdatensätze, die aus den mit dem Influenza-A-Virus infizierten humanen Makrophagenzellen hergestellt wurden, verwendet, um circRNA-Detektion und funktionelle Analysen durchzuführen. Wie im Protokollabschnitt spezifiziert, wurde CIRIquant verwendet, um identifizierte circRNAs zu identifizieren und eine DE-Analyse durchzuführen, wobei die Datensätze der RNA-seq-Bibliothek als Eingabe verwendet wurden. Die verwendeten Referenzdateien basieren auf der neuesten Version des menschlichen Genoms (hg38). Tabelle 4 zeigt ein Beispiel für das endgültige Ergebnis der CIRIquant-Analyse. Die Identifizierung und Filterung von DE-circRNAs aus der CIRIquant-Ausgabe erfolgte durch einfache RStudio-Skripte (Supplementary File 1). CircRNAs werden nur dann als DE klassifiziert, wenn der Wert für die False-Discovery-Rate (FDR) <0,05 und die LogFC-Änderung (LogFC) >|2| beträgt. Tabelle 5 zeigt die Gesamtzahl der nachgewiesenen circRNAs und DE-circRNAs. Insgesamt wurden 35.846 circRNAs nachgewiesen, davon 306 DE. Die in diesem Output detektierten DE-circRNAs werden vollständig hochreguliert (LogFC > 2), wobei keine herunterreguliert wird (LogFC < 2).

Annotation und Charakterisierung von DE-circRNAs
Annotationsstatus von DE-circRNAs
Die identifizierten DE-circRNAs wurden mit einer etablierten circRNA-Datenbank, circBase, abgeglichen. Da die in circBase hinterlegten circRNA-Koordinaten jedoch auf einer früheren humanen Genomversion (hg19) basieren, müssen die circRNA-Koordinaten von circBase für den Abgleich in dieser Studie in die aktuelle Humangenomversion (hg38) konvertiert werden. Zusätzlich muss die Startkoordinate aus der 1-basierten Ausgabe von CIRIquant in 0-basiert konvertiert werden. Die in die hg38-Version konvertierten circRNA-Koordinaten von circBase werden in einem Laufwerksordner in Github (https://github.com/bicciatolab/Circr)³⁷ bereitgestellt. Anschließend wurden die Rstudio-Skripte (Supplementary File 1) verwendet, um den Annotationsstatus von circRNAs in einer neuen Datenrahmenspalte zuzuweisen. Tabelle 6 zeigt ein Beispiel für circRNAs mit dem Annotationsstatus.

Charakterisierung von DE-circRNAs
Dieser Teil wurde vollständig über R-Skripte in der RStudio-Software ausgeführt. R-Skripts vereinfachen die Analyseprozesse, und es sind nur grundlegende Kenntnisse erforderlich.

CircRNA-Typen
In diesem Schritt wurden DE-circRNAs anhand ihrer circRNA-Typen (Antisense, Exonic, Intergenic und Intronic) anhand ihrer genomischen Positionen charakterisiert. Tabelle 7 zeigt die prozentuale Aufschlüsselung der verschiedenen circRNA-Typen, die von den identifizierten DE-circRNAs umfasst werden. Von den insgesamt 306 DE-circRNAs wurden 263 circRNAs (85,95%) als exonische circRNAs identifiziert, der am häufigsten identifizierte circRNA-Typ. Intronische circRNAs sind der am zweithäufigsten identifizierte circRNA-Typ, der 17 DE-circRNAs umfasst, was bis zu 5,56 % der gesamten DE-circRNAs ausmacht. Es folgen intergene circRNAs (16 DE circRNAs ~5,23%) und Antisense circRNAs (10 DE circRNAs ~3,27%).

Anzahl der Gene pro circRNA
CircRNAs, die von CIRIquant identifiziert wurden, können sich über eine Reihe von Genen hinweg überlappen. Bisher konzentrieren sich die meisten Studien auf circRNAs, die sich über ein Gen erstrecken. Daher werden in diesem Protokoll die circRNA-Kandidaten, die mehr als ein Gen umfassen, von der nachgelagerten Analyse ausgeschlossen. Tabelle 8 beschreibt die Anzahl und den prozentualen Anteil der DE-circRNAs, die sich über ein und mehrere Gene erstrecken. In dieser Tabelle werden intergene circRNAs (16 DE-circRNAs) ausgeschlossen, da sie keine Wirtsgene überlappen, während die übrigen circRNA-Typen (290 DE-circRNAs) dieser Analyse unterzogen werden. Von den 290 DE-circRNAs umfasst die Mehrheit der DE-circRNAs (261 circRNAs ~90%) nur ein Gen, während die restlichen 29 circRNAs (~10%) mehr als ein Gen umfassen.

Aufbau des ceRNA-Netzwerks
Ein ceRNA-Netzwerk wird in der Regel gezeichnet, um die circRNA-miRNA-Interaktionen zu visualisieren, nachdem sie vorhergesagt wurden. In Abbildung 3 unten wurde nur eine DE-circRNA als repräsentatives Ergebnis ausgewählt, nämlich die hsa_DE_58 circRNA. Basierend auf Circr-Vorhersagen können hsa_DE_58 bis zu neun verschiedene miRNAs mit einem Schwamm versehen. Diese neun miRNAs werden identifiziert, nachdem sie durch strenge Kriterien gefiltert wurden.

Analyse der funktionalen Anreicherung
GO- und KEGG-Analyse der circRNA-Elterngene
Abbildung 4 unten zeigt ein Blasendiagramm der funktionellen Anreicherung von DE circRNA Elterngenen durch die GO-Analyse. Grundsätzlich zielt die GO-Analyse darauf ab, die biologischen Prozesse, zellulären Standorte und molekularen Funktionen zu entschlüsseln, die bei der untersuchten Erkrankung, in diesem Fall der virusinfizierten Probe, angereichert oder beeinflusst werden. Die Anreicherung wird nur dann als statistisch signifikant betrachtet und im Blasendiagramm dargestellt, wenn der p-Wert < 0,01 liegt. Wie in Abbildung 4 gezeigt, umfassen die drei wichtigsten Anreicherungen für die biologischen Prozesse (BP) die Biogenese des Ribonukleoproteinkomplexes, die Reaktion auf das Virus und die Regulation der Reaktion auf einen biotischen Stimulus, während für die molekularen Funktionen (MF) nur die katalytische Aktivität, die auf RNA wirkt, und die Bindung einzelsträngiger RNA statistisch angereichert sind. Dagegen ist nur der Retromer-Komplex statistisch für die zellulären Komponenten (CC) angereichert.

Abbildung 5 zeigt die KEGG-Anreicherungsanalyse der DE circRNA Elterngene in einem Bubble Plot. Ähnlich wie bei der GO-Anreicherungsanalyse wird die KEGG-Anreicherung nur dann als statistisch signifikant betrachtet und in einem Blasendiagramm dargestellt, wenn der p-Wert < 0,01 beträgt. In diesem Fall wurden nur zwei KEGG-Begriffe angereichert, nämlich der Influenza-A- und der virale Lebenszyklus (HIV-1).

Abbildung 1: Die Pipeline zur Vorhersage und funktionellen Charakterisierung von circRNAs. Die Pipeline zeigt einen einfachen Überblick über die wichtigsten Schritte von Anfang bis Ende, einschließlich der Installation der erforderlichen Softwarepakete, der Vorhersage und Quantifizierung der circRNA-Expression, des Aufbaus des ceRNA-Netzwerks und der Durchführung der funktionellen Anreicherung des circRNA-Elterngens. Bitte klicken Sie hier, um eine größere Version dieser Abbildung zu sehen.

Abbildung 2: Ordnerbaumstruktur für Circr. Diese Ordnerbaumstruktur muss vor dem Ausführen der Circr-Software festgelegt werden, um die erforderlichen Dateien für die Analyse zu erkennen. Bitte klicken Sie hier, um eine größere Version dieser Abbildung zu sehen.

Abbildung 3: ceRNA-Netzwerk, bestehend aus der circRNA-miRNA-Interaktion. Die blaue ovale Form stellt die circRNA dar, während die orangefarbenen Dreiecke die miRNAs darstellen. Die durchgezogenen Linien, die die circRNA mit miRNAs verbinden, beschreiben die potentielle miRNA-Schwammfunktion der hsa_DE_58 circRNA. Bitte klicken Sie hier, um eine größere Version dieser Abbildung zu sehen.

Abbildung 4: Blasendiagramm der GO-Anreicherungsanalyse von DE circRNA-Elterngenen. GeneRatio auf der x-Achse ist die Anzahl der Gene in der Eingabeliste, die dem angegebenen GO-Term zugeordnet sind, und dividiert die Gesamtzahl der Eingabegene. Die Punktgröße im Diagramm wird durch den Zählwert dargestellt, d. h. die Anzahl der Gene in der Eingabeliste, die dem angegebenen GO-Term zugeordnet sind. Je größer die Punkte sind, desto größer ist die Anzahl der Eingabegene, die mit dem Begriff verbunden sind. Außerdem sind die Punkte im Diagramm basierend auf dem p-Wert farbcodiert. Der p-Wert wird berechnet, indem die beobachtete Häufigkeit eines Annotationsterms mit der zufällig erwarteten Häufigkeit verglichen wird. Die einzelnen Terme gelten über einen Cut-off-Wert hinaus als angereichert (p-Wert < 0,01). Der Farbverlauf des p-Wertes von blau bis rot deutet auf eine zunehmende Anreicherung der Terme hin. Bitte klicken Sie hier, um eine größere Version dieser Abbildung zu sehen.

Abbildung 5: KEGG-Anreicherungsanalyse von DE circRNA-Elterngenen. GeneRatio auf der x-Achse ist die Anzahl der Gene in der Eingabeliste, die dem angegebenen KEGG-Term zugeordnet sind, und dividiert die Gesamtzahl der Eingabegene. Die Punktgröße im Diagramm wird durch den Zählwert dargestellt, d. h. die Anzahl der Gene in der Eingabeliste, die dem angegebenen KEGG-Term zugeordnet sind. Je größer die Punkte sind, desto größer ist die Anzahl der Eingabegene, die mit dem Begriff verbunden sind. Außerdem sind die Punkte im Diagramm basierend auf dem p-Wert farbcodiert. Der p-Wert wird berechnet, indem die beobachtete Häufigkeit eines Annotationsterms mit der zufällig erwarteten Häufigkeit verglichen wird. Einzelne Terme gelten als über einen Cut-off-Wert hinaus angereichert (p-Wert < 0,01). Der Farbverlauf des p-Wertes von blau bis rot deutet auf eine zunehmende Anreicherung der Terme hin. Bitte klicken Sie hier, um eine größere Version dieser Abbildung zu sehen.

Tabelle 1: Die .lst-Dateivorbereitung von CIRIquant. Die Zielpfade der Kontroll- und behandelten Proben aus der CIRIquant-Ausgabe werden in eine Textdatei geschrieben, um die Expressionen der circRNA zwischen den beiden Probentypen zu vergleichen.

Tabelle 2: Beispiel für eine BED-Datei für Circr. Sechs Spalten (Chr, Start, Ende, Name, Gen und Strang), die mit den circRNAs verknüpft sind, sind erforderlich, um die BED-Datei zu generieren.

Tabelle 3: Cytoscape-Eingabedatei. Vier Spalten (circRNA_name, Type, miRNA_name und Type) müssen in eine Textdatei geschrieben werden.

Tabelle 4: Teil der finalen Ausgabedatei (.csv) von CIRIquant. CIRIquant liefert Informationen wie LogFC, Log-Counts per Million (LogCPM), logistische Regression (LR), p-Wert, Differentialausdruck und FDR.

Tabelle 5: Eine Zusammenfassung der Gesamtzahl der identifizierten und differentiell exprimierten (DE) circRNAs. Insgesamt wurden 35.846 circRNAs nachgewiesen, davon 306 DE-circRNAs. Alle 306 DE-circRNAs sind in behandelten Proben im Vergleich zu Kontrollproben hochreguliert (wobei keine herunterreguliert wird).

Tabelle 6: Tabelle der benutzerdefinierten circRNA-Namen mit Annotationsstatus. CircRNAs werden in einer Datenbank bekannter hinterlegter circRNAs (circBase) abgefragt. Wenn die circRNA in der Datenbank vorhanden ist, wird sie als annotiert markiert, während das Fehlen der circRNA als nicht annotiert markiert wird.

Tabelle 7: Identifizierte circRNA-Typen. CircRNAs können basierend auf ihrer Sequenzregion weiter in verschiedene Arten von circRNAs kategorisiert werden, nämlich exonisch, intronisch, antisense und intergen.

Tabelle 8: Prozentualer Anteil der circRNAs mit der unterschiedlichen Anzahl von Genen überspannt. CircRNAs werden üblicherweise von Exons eines Gens kodiert, aber auch circRNAs, die mehr als ein Gen umfassen, können von CIRIquant nachgewiesen werden.

Ergänzende Datei 1: Skripte, die im Protokoll verwendet werden. Bitte klicken Sie hier, um diese Datei herunterzuladen.

Um den Nutzen dieses Protokolls zu veranschaulichen, wurde RNA-seq aus mit dem Influenza-A-Virus infizierten humanen Makrophagenzellen als Beispiel verwendet. CircRNAs, die als potentielle miRNA-Schwämme in Wirt-Pathogen-Interaktionen fungieren, und ihre GO- und KEGG-funktionelle Anreicherung innerhalb eines Wirts wurden untersucht. Obwohl es eine Vielzahl von circRNA-Tools gibt, die online verfügbar sind, ist jedes von ihnen ein eigenständiges Paket, das nicht miteinander interagiert. Hier stellen wir einige der Werkzeuge zusammen, die für die circRNA-Vorhersage und -Quantifizierung, die funktionelle Anreicherung von circRNAs, die Vorhersage der circRNA-miRNA-Interaktion und den Aufbau von ceRNA-Netzwerken erforderlich sind. Dieses optimierte Protokoll ist zeitsparend und kann auf klinische Proben angewendet werden, um circRNA-Kandidaten mit diagnostischen und prognostischen Werten nachzuweisen.

Im Wesentlichen haben wir CIRIquant³¹ verwendet, ein circRNA-Quantifizierungstool, das mit CIRI2 vorverpackt ist und die DE-Analyse von circRNAs nachweisen und durchführen kann. DE circRNAs werden auf der Grundlage eines Cut-off-Wertes von LogFC > |2| und FDR < 0,05, was dazu beiträgt, potenzielle Fehlalarme in nachgelagerten Analysen zu eliminieren. Die Charakterisierung von DE-circRNAs in Bezug auf den Annotationsstatus, die circRNA-Typen und die Anzahl der gestreckten Gene hilft bei der Kategorisierung und weiteren Filterung von circRNA-Kandidaten. Anschließend wird Circr³⁷, ein circRNA-miRNA-Vorhersagewerkzeug, verwendet, um potenzielle miRNA-Schwammkandidaten vorherzusagen. Nach der Vorhersage potenzieller miRNAs als Ziele von circRNAs wird ein ceRNA-Netzwerk gezeichnet. Schließlich wird auf der Grundlage der elterlichen Gene von circRNAs das R clusterProfiler-Paket³⁹ für die funktionelle Annotation über die GO- und KEGG-Signalweg-Anreicherungsanalyse verwendet. Die Ergebnisse von GO und KEGG könnten dazu beitragen, die biologischen Mechanismen zu entschlüsseln, die von circRNAs beeinflusst werden.

Bis heute wurden verschiedene circRNA-Vorhersagetools entwickelt, darunter CIRI2 43, CIRCexplorer2 44, find_circ 45, MapSplice 46 und UROBORUS ⁴⁷. In einer von Hansen et al. durchgeführten Studie wird berichtet, dass CIRI2 eine hohe Gesamtleistung aufweist. Es gehört zu den wenigen circRNA-Nachweisinstrumenten, die im Hinblick auf die De-novo-Vorhersage und die Reduzierung der Identifizierung falsch positiver Ergebnisse gut funktionieren^{können 48}. In dieser Studie wurde daher CIRIquant verwendet, das CIRI2 für die circRNA-Detektion und -Quantifizierung verwendet. CIRIquant wurde verwendet, um die BSJ-Reads (Back Splice Junction) zu zählen, und die Zähldaten wurden auf die Reads normalisiert, die auf verwandte lineare RNAs abgebildet wurden, die von denselben Genloci transkribiert wurden. Dies ermöglicht die Quantifizierung von circRNAs in einer Probe. Um die differentielle Expression von circRNAs unter experimentellen Bedingungen zu bestimmen, implementierte CIRIquant ein verallgemeinertes lineares Modell in edgeR⁴⁹ für die DE-Analyse, und der exakte Rate-Ratio-Test wurde als statistischer Test verwendet, um die Signifikanz des Unterschieds im circRNA-Übergangsverhältnis zu bestimmen. Obwohl andere circRNA-Quantifizierungsinstrumente wie CIRCexplorer3-CLEAR⁵⁰ verwendet werden können, um das Expressionsniveau von circRNAs zu quantifizieren, erlaubt dieses Tool nur die circRNA-Quantifizierung in einer Probe, da es die BSJ-Lesevorgänge in einer Probe zählt und die Zähldaten mit den verwandten linearen RNA-Zählungen derselben Probe vergleicht. CIRCexplorer3-CLEAR kann circRNA-Expressionen unter experimentellen Bedingungen nicht vergleichen. Darüber hinaus ist in CIRCexplorer3-CLEAR kein statistisches Analysewerkzeug implementiert, um die quantifizierte Expressionsebene zu unterstützen. Obwohl das standardmäßige circRNA-Vorhersagewerkzeug, das in CIRIquant implementiert ist, CIRI2 ist, können die Vorhersageergebnisse von anderen Werkzeugen wie find_circ und CIRCexplorer2 auch für die Quantifizierung und DE-Analyse³¹ verwendet werden. In diesem Protokoll wurde nur ein circRNA-Vorhersagewerkzeug (CIRI2) für die Vorhersage verwendet, das immer noch falsch-positive circRNA-Kandidaten liefern könnte. Um falsch positive Ergebnisse zu reduzieren, kann man andere circRNA-Vorhersageinstrumente für die Analyse kombinieren und gemeinsame circRNAs auswählen, die unter den verschiedenen circRNA-Vorhersagewerkzeugen entdeckt wurden^48,51. Um die circRNA-Detektion weiter zu verbessern, ist es ideal, RNA-Sequenzierungsdatensätze zu verwenden, die sowohl rRNA-depletiert als auch einer RNase R-Vorbehandlung unterzogen werden.

Je nach Ziel der Studie können de novo und annotierte DE circRNAs auf Basis der circBase-Datenbank separat identifiziert werden⁵². CirrcRNAs, die mehr als ein Gen umfassen, müssen jedoch häufig manuell auf UCSC oder einem anderen Genombrowser untersucht werden, um die Authentizität von circRNAs zu bestimmen und falsch positive Ergebnisse zu eliminieren. Nichtsdestotrotz wurden kürzlich auch circRNAs beschrieben, die mehr als ein Gen umfassen, wie z. B. circRNAs, die von Fusionsgenen abgeleitet sind^53,54.

Circr kombiniert drei verschiedene miRNA-mRNA-Vorhersagealgorithmen, nämlich TargetScan⁵⁵, miRanda 56 und RNAhybrid⁵⁷, um die circRNA-miRNA-Bindungsstellen vorherzusagen. Darüber hinaus bezieht der Algorithmus auch Informationen über AGO-Peaks und zuvor validierte Wechselwirkungen in die circRNA-miRNA-Analyse ein. Hier wurden strenge Filterkriterien angewendet, um eine zuverlässigere circRNA-miRNA-Vorhersage zu erhalten und so die Zahl der falsch positiven Ergebnisse weiter zu reduzieren. Die Stringenz dieses Filterschritts kann jedoch je nach Benutzerpräferenz höher oder niedriger eingestellt werden.

ClusterProfiler ist ein gut dokumentiertes R-Paket, das Gensätze in verschiedenen Organismen funktionell annotieren kann. Neben den Funktionen innerhalb des R clusterProfiler-Pakets, die in diesem Protokoll erwähnt werden (enrichGO und enrichKEGG), die die Überrepräsentationsanalyse verwenden, gibt es auch andere Funktionen wie gseGO und gseKEGG, die verwendet werden können. Wenn clusterProfiler keine geeignete Wahl für den Workflow ist, gibt es auch andere Tools und Pakete wie den "AllEnricher"⁵⁸ oder die websitebasierten Tools wie "Metascape"⁵⁹, die eine Reihe von Genen funktional annotieren können. Obwohl die oben beschriebene Pipeline bei der Vorhersage potenzieller circRNAs und ihrer funktionellen Annotationen hilft, ist eine Nasslaborverifizierung erforderlich, um solide Beweise zu liefern.

Die Autoren haben nichts zu verraten.

Der Autor dankt Tan Ke En und Dr. Cameron Bracken für die kritische Durchsicht dieses Manuskripts. Diese Arbeit wurde durch Zuschüsse aus dem Fundamental Research Grant Scheme (FRGS/1/2020/SKK0/UM/02/15) und dem University of Malaya High Impact Research Grant (UM. C/625/1/HIR/MOE/CHAN/02/07).