De conceptuele schets van GREIN wordt getoond in Fig. 1. Individuele RNA-seq datasets worden verwerkt door de GREP2 pipeline en lokaal opgeslagen als R Expression Sets. Gebruikers hebben toegang tot voorbewerkte datasets en kunnen deze analyseren via de GREIN grafische gebruikersinterface (GUI) of kunnen datasets die nog niet verwerkt zijn, indienen voor verwerking. GUI-gestuurde workflows vergemakkelijken het onderzoek en de visualisatie van data, statistische analyse, constructie van transcriptionele handtekeningen, en systeembiologische interpretatie van differentieel uitgedrukte (DE) genen. Zowel GREIN als de back-end pijplijn (GREP2) zijn geschreven in R en uitgebracht als respectievelijk Docker container en R package. Grafische gebruikersinterfaces voor GREIN zijn geïmplementeerd in Shiny16, een web framework voor het bouwen van dynamische web applicaties in R. De web instantie op https://shiny.ilincs.org/grein wordt ingezet via een robuuste Docker zwerm van load-balanced Shiny servers. De volledige GREIN-infrastructuur, inclusief de verwerkingspijplijn, wordt ingezet via Docker-containers.
Schematische workflow van GREP2, webinterface en outputs van GREIN. GEO-datasets worden systematisch verwerkt met de GREP2-pijplijn en opgeslagen in de back-end datasetbibliotheek. GUI-gestuurde GREIN-workflows vergemakkelijken uitgebreide analyse en visualisatie van verwerkte datasets.
Gebruikersvriendelijke GUI-gestuurde workflows in GREIN vergemakkelijken typische hergebruikscenario’s voor RNA-seq-gegevens, zoals onderzoek van kwaliteitscontrolemaatregelen en visualisatie van expressiepatronen in de hele dataset, steekproefgrootte en power-analyse ter informatie van experimenteel ontwerp van toekomstige studies, statistische differentiële genexpressie, verrijking van genenlijsten en netwerkanalyse. Naast de standaard twee-groepen vergelijking, ondersteunt de differentiële genexpressie analyse module ook de aanpassing van een gegeneraliseerd lineair model dat rekening houdt met covariaten of batch-effecten. De interactieve visualisatie en exploratie tools geïmplementeerd omvatten clusteranalyse, interactieve heatmaps, principale componenten analyse (PCA), t-verdeelde stochastische buur embedding (t-SNE), enz. (Aanvullende tabel S1). Gebruiker kan ook zoeken naar ontologische annotaties van menselijke RNA-seq monsters en datasets die door de MetaSRA project 13. Elke verwerkte menselijke RNA-seq monster is gelabeld met MetaSRA in kaart brengen van biomedische ontologieën, waaronder Disease Ontology, Cell Ontology, Experimental Factor Ontology, Cellosaurus, en Uberon. Biologische interpretatie van differentiële genexpressies wordt geholpen door directe links naar andere online tools voor het uitvoeren van typische post-hoc analyses, zoals de genenlijst en pathway verrijkingsanalyse en de netwerkanalyse van differentieel tot expressie komende (DE) genen. De verbinding met deze analytische webdiensten wordt uitgevoerd door het indienen van de differentiële genexpressie handtekening (dat wil zeggen, de lijst van gemiddelde veranderingen in genexpressie en bijbehorende p-waarden voor alle / up / down gereguleerde genen geanalyseerd) naar iLINCS17 (Integrative LINCS). iLINCS biedt ook de handtekeningen connectiviteit analyse voor de onlangs vrijgegeven Connectivity Map L1000 handtekeningen18. Gedetailleerde stap-voor-stap instructies over GREIN analyse workflows worden gegeven in het aanvullend materiaal en ‘Help’ sectie in GREIN.
Key functionaliteiten
Zoeken of indienen voor verwerking
Gebruiker kan ofwel zoeken naar een reeds verwerkte GEO-dataset in het vak ‘Zoeken naar GEO-serie (GSE) toetreding’ of een dataset indienen voor verwerking als de dataset nog niet is verwerkt (Supplementary Fig. S2). Op dit punt in de tijd, de overgrote meerderheid van GEO menselijke, muis, en rat RNA-seq datasets zijn voorbewerkt en de gebruiker-submission van GEO datasets voor de verwerking zal worden vereist slechts af en toe. Gebruiker kan de verwerking status van de gevraagde dataset in de ‘Processing console’ tab (Supplementary Fig. S3). Andere zoekopties omvatten het zoeken op trefwoorden via metadata van de datasets en het zoeken van monsters via biomedische ontologieën via MetaSRA ontologische annotaties.
Explore dataset
GREIN biedt toegang tot zowel ruwe als genormaliseerde (tellingen per miljoen en transcript per miljoen) gen- en transcriptniveaugegevens. GREIN wordt geleverd met verschillende interactieve en aanpasbare tools om expressiepatronen te visualiseren, zoals interactieve heatmaps van geclusterde genen en monsters, dichtheid plots voor alle of een subset van monsters, tussen en binnen groep variabiliteit analyse door middel van 2D en 3D dimensionaliteit reductie analyses en visualisaties zoals PCA en t-SNE (Fig. 2). De gebruiker kan ook het expressieprofiel van elk gen afzonderlijk visualiseren (Supplementary Fig. S6).
Exploratieve analyseplots in GREIN. (A) Correlatie heatmap toont een hogere correlatie binnen cellijnen en een lage correlatie tussen cellijnen. Over het algemeen duiden hoge correlaties binnen elke cellijn op een hoge kwaliteit van de transcriptionele profielen. (B) Hiërarchische clustering op basis van Pearson correlatie van de top 500 meest variabele genen op basis van de mediaan absolute afwijking als de variabiliteit maatregel. De gegevens zijn genormaliseerd en gecentreerd naar het gemiddelde. (C) Drie-dimensionale principale componenten analyse plot van de cellijnen. (D) Twee-dimensionale t-SNE plot van behandeling voorwaarde en cellijn toont duidelijke scheiding van de cellijnen, en vervolgens de RNA fracties met vermelding van twee dominante bronnen van de variabiliteit tussen RNA-seq profielen.
Quality control
De kwaliteit van de RNA-seq gegevens in openbare repositories blijft een groot probleem. In een recente studie van Deelen et al.19 moest meer dan de helft van de 65.000 verwerkte openbare RNA-seq-monsters worden verwijderd als gevolg van QC-problemen. In plaats van het verwijderen van monsters, GREIN biedt een uitgebreide kwaliteitscontrole (QC) verslag van de ruwe sequentie data en sequentie mapping voor elk monster (Supplementary Fig. S7), en stelt de gebruiker in staat om een beslissing te nemen over welke monsters moeten worden uitgesloten van downstream analyses.
Statistische power analyse
De power analyse module in GREIN vergemakkelijkt de berekening en visualisatie van de statistische kracht van het opsporen van differentieel tot expressie komende genen in toekomstige studies met behulp van vergelijkbare biologische monsters. Het schatten van de juiste steekproefgrootte voor toekomstige studies met vergelijkbare biologische monsters is vaak de belangrijkste motiverende factor in de heranalyse van RNA-seq gegevens. Power analyse vergemakkelijkt ook de post-hoc analyse van vals negatieve percentages in de huidige dataset. Het gebrek aan statistisch vermogen en verschillen in statistisch vermogen tussen genen kunnen vals-negatieve resultaten opleveren die tot verkeerde conclusies leiden20. De ‘Power curve’ segment biedt macht schattingen voor verschillende aantal monsters op basis van een enkel gen (Fig. 3A). De gebruiker kan de standaardwaarden van de parameters wijzigen. De plot ‘Detecteerbaarheid van genen’ visualiseert de vermogensschatting van elk van de genen op basis van de geselecteerde groepen en gen-wise dispersie (Fig. 3B). De gemiddelde dekking van de genen wordt uitgezet tegen hun biologische variabiliteit en wordt weergegeven in twee sets op basis van hun detecteerbaarheidsstatus (vermogen ≥0,8 en vermogen < 0,8).
Machtsanalyse voor het beoordelen van transcriptionele veranderingen in niet-maligne MCF10A-cellijn. (A) Schattingen van het vermogen op basis van één gen voor een verschillend aantal monsters in elke groep met een minimale vouwverandering van 2 en statistische significantie α = 0,01. (B) Gen-wise detecteerbaarheid op het log2CPM-BCOV vlak met FDR ≤0.1 en twee monsters in elke groep.
Differentiële genexpressie
Het creëren en interpreteren van differentiële genexpressie handtekening is een typisch analyse scenario in RNA-seq experimenten. Met GREIN kan de gebruiker een handtekening creëren door genexpressie tussen twee groepen monsters te vergelijken, met of zonder aanpassingen voor experimentele covariaten of batch-effecten. GREIN kan complexe experimentele ontwerpen aan door de flexibiliteit te bieden van het herschikken van groepen en subgroepen of het selecteren van specifieke monsters. Differentiële expressie handtekening kan worden gevisualiseerd via interactieve grafieken die heatmap van top differentieel gedereguleerde genen (Supplementary Fig. S15) gerangschikt naar false discovery rate (FDR), log fold change vs. log gemiddelde expressie (MA) plot (Supplementary Fig. S16), en gen detectability plot (Supplementary Fig. S17) omvatten. Differentiële expressie handtekening, met of zonder boekhouding voor potentieel vals-negatieve resultaten, kan direct worden geëxporteerd naar iLINCS voor verrijking en connectiviteit analyse.
Use case: Analyse van transcriptionele en translationele regulatie van hypoxie in niet-maligne borst epitheliale en triple-negatieve borstkanker cellijnen
We demonstreren het gebruik van GREIN door het opnieuw analyseren van een onlangs gepubliceerde GEO RNA-seq gegevens (GSE104193). Sesé et al.21 onderzochten de transcriptionele en translationele regulatie van hormoon-refractaire triple-negatieve borstkanker (TNBC) subtype onder een combinatie van hypoxie en mTOR (mechanistic target of rapamycin) remmer behandeling. In het bijzonder analyseerden de auteurs de expressieprofielen van TNBC (MDA-MB-231) en niet-maligne borst epitheliale (MCF10A) cellen blootgesteld aan normoxische (21% O2) en hypoxische (0,5% O2) omstandigheden en/of behandeld met een mTORC1 en -2 inhibitor PP242. Elk van de monsters werd gesequenced voor totaal (T) en polysoom-gebonden (P) mRNA. De dataset bevat 32 monsters, die twee biologische replicaten voor elke combinatie van cellijn, zuurstof niveau, behandeling status, en mRNA fractie.
Exploratieve analyse van de verwerkte dataset in GREIN (Fig. 2) blijkt dat de sterkste bron van variatie in tussen de monsters komt van verschillen tussen de twee cellijnen. Dit wordt bevestigd door de correlatie analyse van volledige expressie profielen (Fig. 2A), de hiërarchische clustering van top 500 zeer variabele genen op basis van mediane absolute afwijking (Fig. 2B), 3D PCA plot van de monsters (Fig. 2C), en de 2D t-SNE plot (Fig. 2D). Bovendien hoge correlaties tussen expressie profielen voor dezelfde cellijn (Fig. 2A) wijst op een goede signaal-ruis in de genexpressie metingen. De extra substructuur van gegevens aangegeven door de 2D t-SNE plot is onderzocht door het schilderen van monsters op basis van verschillende attributen (Supplementary Fig. S1). Deze analyse toonde aan dat scheidingen binnen elke cellijn worden geïnduceerd door verschillende mRNA fracties en vervolgens verschillen tussen experimentele condities.
Next, gebruikten we GREIN om statistische power analyse uit te voeren op basis van het patroon van biologische variabiliteit waargenomen in deze dataset. We overwogen transcriptionele profielen van elke cellijn blootgesteld aan hypoxie en behandeld met of zonder PP242, wat leidt tot vier vergelijkingen. Uitgaande van een expressieverschil van ten minste twee maal tussen de groepen, bij de statistische significantie van α = 0,01, en met slechts twee herhalingen in elke groep, is de statistische power van een gen om te worden gedetecteerd als differentieel tot expressie gebracht lager dan 0,55 in alle vergelijkingen (Tabel 2). Onze analyse geeft aan dat men vier herhalingen per groep nodig zou hebben om 80% vermogen te bereiken bij het detecteren van tweevoudige verandering in expressie (Tabel 2 en Fig. 3A). In een typische RNA-seq experiment, een sequencing diepte van 20-30 miljoen is voldoende om genexpressie kwantificeren voor bijna alle genen 4,22 die ook duidelijk is in deze dataset. We hebben ook geëvalueerd statistisch vermogen van elk gen te worden gedetecteerd als differentieel tot expressie uit de ‘Detecteerbaarheid van genen’ plot. De gemiddelde log van counts per million (CPM) waarden van de genen werden uitgezet tegen de gen-wise biologische variatiecoëfficiënt (BCOV) en de power werd berekend voor de corresponderende genen (Fig. 3B). Een gecontroleerde false discovery rate van 0,05 en een verwacht percentage ware positieven van 10% werden gebruikt om de statistische significantie te schatten. We definiëren een gen als detecteerbaar als differentieel tot expressie komend in hypoxische conditie als de power 0.8 of hoger is. Zoals verwacht bestaat er een omgekeerd evenredige relatie tussen BCOV en vermogen (Fig. 3B). Ook neemt het vermogen om differentiële expressie van een gen te detecteren toe met een hogere log CPM of effect size.
Eén van de doelen van de studie was om transcriptionele veranderingen in hypoxische en normoxische condities met en zonder PP242 behandeling in zowel MCF10A als MDA-MB-231 cellijnen te analyseren. We creëerden transcriptionele handtekeningen van hypoxie en hypoxie + PP242 in totaal mRNA door differentiële expressie analyse tussen hypoxie en hypoxie + PP242 monsters respectievelijk ten opzichte van de controlemonsters, terwijl we corrigeerden voor batch effect door ‘repliceren’ te behandelen als een covariaat, voor elke cellijn afzonderlijk. Wij vonden een hoger aantal genen differentieel tot expressie (DE) in MCF10A cellijnen vergeleken met MDA-MB-231 in zowel hypoxie als hypoxie + PP242 (Fig. 4A), wat erop wijst dat de tumor cellijn misschien beter is toegerust om met hypoxie om te gaan. Deze analyse toonde ook aan dat de meeste niet-differentieel geëxpresseerde genen ook niet detecteerbaar zijn, wat erop wijst dat zij vals-negatieve resultaten kunnen vertegenwoordigen. Dit is in overeenstemming met de power-analyse waaruit blijkt dat 4 monsters per groep nodig zouden zijn om differentieel tot expressie komende genen consistent te identificeren met een gemiddelde BCOV. Om lager tot expressie komende genen te identificeren zou een nog grotere steekproefgrootte nodig zijn.
Differentiële expressie en detecteerbaarheid van de genen. (A) Het aantal genen (log10 schaal) niet differentieel tot expressie en niet detecteerbaar (NDE&NDT), differentieel tot expressie (DE), en niet differentieel tot expressie, maar detecteerbaar (NDE&DT) in de vergelijkingen met normoxie voor de totale mRNA-fractie. We noemen een gen detecteerbaar (DT) als het vermogen ≥0,8 en differentieerbaar als FDR < 0,05. (B) De gen-detecteerbaarheidsplot voor de eerste vergelijking (MCF10A en hypoxie) die de bovengenoemde lijst van genen samen met hun respectieve vouwveranderingen (FC) visualiseert.
Om differentieel tot expressie komende genen te interpreteren in termen van beïnvloede biologische routes, hebben we de differentiële genexpressiehandtekeningen van hypoxie voorgelegd aan online verrijkingshulpmiddelen (DAVID23, ToppGene24, Enrichr25, en Reactome26) via iLINCS. De ingediende handtekeningen omvatten een gecombineerde lijst van DE en NDE&DT genen die waarschijnlijk echte positieve en echte negatieve genen vertegenwoordigen. Genen werden geselecteerd op basis van een cutoff van 0,7 en 0,01 voor statistische power en FDR respectievelijk. Figuur 5 illustreert de verrijkingsresultaten verkregen met ToppGene voor de MCF10 hypoxie signatuur. Significant verrijkt (FDR < 0.05) top 10 gen ontologie (GO) categorieën uit ToppGene en DAVID functionele annotatie-tool omvatten reactie op hypoxie, reactie op verlaagde zuurstofniveaus, angiogenese, regulering van celproliferatie, oxidatie-reductie proces, en reactie op abiotische stimulus die gemeenschappelijk zijn in beide cellijnen (Supplementary Table S2 en Supplementary Table S3). De meeste van deze categorieën zijn consistent met de oorspronkelijke studie. Bovendien identificeerde ToppGene suite hypoxie geïnduceerde factor (HIF-1-alpha) transcriptiefactor netwerk dat werd geactiveerd in beide cellijnen (Supplementary Table S4 en Supplementary Table S5).
Snapshot van enkele van de significante pathway en gen ontologie (GO) categorieën uit ToppGene via iLINCS. Deze categorieën zijn gevonden in de vergelijking tussen hypoxie en normoxie in MCF10A cellijn met behulp van een gecombineerde lijst van DE en NDE&DT genen. De rode verticale lijn is de geselecteerde cutoff van 0,05.
Ten slotte gebruikten we GREIN-verbinding met iLINCS om de geüploade handtekening te “verbinden” met LINCS27-consensus (CGS) gen knockdown handtekeningen18. We vonden 3.727 LINCS consensus gen knockdown handtekeningen die significant (pValue < 0,05) verbonden waren met onze geüploade handtekening. De doelgenen van de top 100 aangesloten handtekeningen werden geselecteerd voor verdere verrijkingsanalyse. We vonden cellulaire respons op hypoxie en regulering van Hypoxia-inducible Factor (HIF) door zuurstof in de lijst van top 10 geactiveerde paden in beide cellijnen (Supplementary Table S6 en Supplementary Table S7). Hoewel deze analyse vergelijkbare verrijkte functionele categorieën oplevert als de initiële verrijkingsanalyse, vult ze de oorspronkelijke analyse aan door verschillende doelgenen te impliceren die niet differentieel tot expressie komen, hoewel ze voldoende tot expressie komen om volgens onze power-analyses detecteerbaar te zijn. Door deze twee resultaten aan elkaar te koppelen, worden deze genen geïmpliceerd als potentiële regulatoren op een hoger niveau van de respons op hypoxie.