O esboço conceptual de GREIN é mostrado na Fig. 1. Os conjuntos de dados individuais do RNA-seq são processados pelo gasoduto GREP2 e armazenados localmente como Conjuntos de Expressão R. O usuário pode acessar e analisar conjuntos de dados pré-processados através da interface gráfica GREIN (GUI) ou submeter para processamento conjuntos de dados que ainda não tenham sido processados. Os fluxos de trabalho guiados pela GUI facilitam o exame e visualização de dados, análise estatística, construção de assinaturas transcripcionais e interpretação da biologia de sistemas de genes diferencialmente expressos (DE). Tanto o GREIN quanto o back-end pipeline (GREP2) são escritos em R e liberados como Docker container e pacote R, respectivamente. As interfaces gráficas de usuário para o GREIN são implementadas em Shiny16, uma estrutura web para construir aplicações web dinâmicas em R. A instância web em https://shiny.ilincs.org/grein é implementada através de um robusto enxame Docker de servidores Shiny balanceados em carga. A infra-estrutura completa do GREIN, incluindo o pipeline de processamento é implantada via Docker containers.
Schematic workflow of GREP2, interface web e saídas do GREIN. Os conjuntos de dados GEO são sistematicamente processados usando o gasoduto GREP2 e armazenados na biblioteca do conjunto de dados back-end. Os fluxos de trabalho guiados por GUI no GREIN facilitam a análise e visualização abrangente dos conjuntos de dados processados.
Os fluxos de trabalho guiados por GUI no GREIN facilitam cenários típicos de reutilização de dados RNA-seq, tais como exame de medidas de controle de qualidade e visualização de padrões de expressão em todo o conjunto de dados, tamanho da amostra e análise de potência com o propósito de informar o desenho experimental de estudos futuros, expressão diferencial de genes estatísticos, enriquecimento da lista de genes e análise de rede. Além da comparação padrão de dois grupos, o módulo de análise de expressão gênica diferencial também suporta o ajuste de um modelo linear generalizado que contabiliza covariates ou efeitos de lote. As ferramentas de visualização e exploração interativa implementadas incluem análise de cluster, heatmaps interativos, análise de componentes principais (PCA), incorporação de vizinhos estocásticos distribuídos em t (t-SNE), etc. (Tabela suplementar S1). O usuário também pode procurar por anotações ontológicas de amostras e conjuntos de dados de RNA-seq humanos fornecidos pelo projeto MetaSRA13. Cada amostra processada de RNA-seq humano é rotulada com MetaSRA mapeamento de ontologias biomédicas incluindo Ontologia da Doença, Ontologia Celular, Ontologia do Fator Experimental, Celosaurus e Uberon. A interpretação biológica de expressões diferenciais de genes é auxiliada por links diretos com outras ferramentas online para realizar análises típicas pós-hoc, como a lista de genes e a análise de enriquecimento do caminho e a análise em rede de genes expressos diferencialmente (DE). A conexão a estes serviços web analíticos é implementada através da submissão da assinatura da expressão diferencial de genes (ou seja, a lista de mudanças médias na expressão gênica e os valores p associados para todos os genes regulados (up/down) analisados) ao iLINCS17 (LINCS integrativo). O iLINCS também fornece a análise de conectividade de assinaturas para assinaturas L1000 do Mapa de Conectividade recentemente lançado18. Instruções detalhadas passo a passo sobre os fluxos de trabalho da análise GREIN são fornecidas na seção Material Suplementar e ‘Ajuda’ em GREIN.
Funcionalidades-chave
Procurar ou submeter para processamento
O usuário pode procurar por um conjunto de dados GEO já processado na caixa ‘Search for GEO series (GSE) accession’ ou submeter um conjunto de dados para processamento se o conjunto de dados ainda não estiver processado (Fig. Suplementar S2). Neste momento, a grande maioria dos conjuntos de dados GEO humanos, mouse e ratazanas RNA-seq foram pré-processados e o envio pelo usuário de conjuntos de dados GEO para processamento só será necessário ocasionalmente. O usuário pode verificar o status do processamento do conjunto de dados solicitado na guia ‘Console de processamento’ (Fig. S3 suplementar). Outras opções de pesquisa incluem pesquisa por palavras-chave através de metadados dos conjuntos de dados e pesquisa de amostras através de ontologias biomédicas através de anotações ontológicas MetaSRA.
Explorar conjunto de dados
GREIN permite o acesso tanto a dados brutos como a dados normalizados (contagens por milhão e transcrição por milhão) do gene e do nível de transcrição. GREIN vem com várias ferramentas interativas e customizáveis para visualizar padrões de expressão como heatmaps interativos de genes e amostras agrupadas, gráficos de densidade para todas ou um subconjunto de amostras, entre e dentro do grupo análise de variabilidade através de análises e visualizações de redução de dimensionalidade 2D e 3D como PCA e t-SNE (Fig. 2). O usuário também pode visualizar o perfil de expressão de cada gene separadamente (Figura Suplementar S6).
Gráficos de análise explorativa em GREIN. (A) O heatmap de correlação mostra uma correlação mais alta dentro das linhas celulares e uma correlação baixa entre as linhas celulares. Geralmente correlações altas dentro de cada linha celular indicam uma alta qualidade de perfis transcripcionais. (B) Agrupamento hierárquico baseado na correlação de Pearson dos 500 genes mais variáveis com base no desvio absoluto mediano como medida de variabilidade. Os dados são normalizados e centrados na média. (C) Gráfico tridimensional de análise dos componentes principais das linhas celulares. (D) O gráfico bidimensional t-SNE da condição de tratamento e linha celular mostra clara separação das linhas celulares, e então as frações de RNA indicando duas fontes dominantes da variabilidade entre os perfis RNA-seq.
Controle de qualidade
A qualidade dos dados do RNA-seq em repositórios públicos continua a ser um grande problema. Em um estudo recente de Deelen et al.19, mais da metade das 65.000 amostras de RNA-seq públicas processadas tiveram que ser removidas devido a problemas de QC. Ao invés de remover amostras, o GREIN fornece um relatório completo de controle de qualidade (QC) de dados de seqüência bruta e mapeamento de seqüência para cada amostra (Figura Complementar S7), e permite ao usuário tomar uma decisão sobre quais amostras devem ser excluídas das análises posteriores.
Análise de potência estatística
O módulo de análise de potência no GREIN facilita o cálculo e visualização do poder estatístico de detecção de genes expressos diferentemente em estudos futuros utilizando amostras biológicas similares. Estimar o tamanho apropriado da amostra para estudos futuros com amostras biológicas similares é muitas vezes o principal fator motivador na re-análise de dados RNA-seq. A análise de poder também facilita a análise post-hoc de taxas de falsos negativos no conjunto de dados atual. A falta de poder estatístico e as diferenças no poder estatístico entre genes podem produzir resultados falsos negativos levando a conclusões erradas20. O segmento ‘Power curve’ fornece estimativas de potência para diferentes números de amostras com base em um único gene (Fig. 3A). O usuário pode modificar os valores padrão dos parâmetros. O gráfico ‘Detectabilidade dos genes’ visualiza a estimativa de potência de cada um dos genes com base nos grupos seleccionados e na dispersão no sentido do género (Fig. 3B). A cobertura média dos genes é plotada contra sua variabilidade biológica e é exibida em dois conjuntos baseados em seu estado de detectabilidade (potência ≥0.8 e potência < 0.8).
Análise de potência para avaliar mudanças transcripcionais na linha de células MCF10A não malignas. (A) Estimativas de potência baseadas em um único gene para diferentes números de amostras em cada grupo com uma mudança mínima de dobra de 2 e significância estatística α = 0,01. (B) Detectabilidade no plano log2CPM-BCOV com FDR ≤0.1 e duas amostras em cada grupo.
Expressão diferencial de genes
Criar e interpretar a assinatura de expressão diferencial de genes é um cenário de análise típico em experimentos RNA-seq. Com GREIN, o usuário pode criar uma assinatura comparando a expressão gênica entre dois grupos de amostras com ou sem ajustes para covariates experimentais ou efeitos de lote. O GREIN pode lidar com desenhos experimentais complexos fornecendo a flexibilidade de rearranjar grupos e subgrupos ou selecionar amostras específicas. A assinatura da expressão diferencial pode ser visualizada através de gráficos interativos que incluem o mapa térmico dos principais genes diferencialmente desregulados (Fig. Complementar S15) classificados por taxa de falsas descobertas (FDR), gráfico de log fold change vs. log average expression (MA) (Fig. Complementar S16), e gráfico de detecção de genes (Fig. Complementar S17). A assinatura de expressão diferencial, com ou sem contabilização de resultados potencialmente falsos negativos, pode ser exportada diretamente para o iLINCS para análise de enriquecimento e conectividade.
Caso de utilização: Análise da regulação transcricional e translacional da hipoxia em linhas de células não malignas do epitélio mamário e tri-negativas do cancro da mama
Demonstramos o uso do GREIN, reanalisando um dado GEO RNA-seq recentemente publicado (GSE104193). Sesé et al.21 examinaram a regulação transcricional e translacional do subtipo de câncer de mama tríplice negativo com hormônio refratário (TNBC) sob uma combinação de hipóxia e mTOR (alvo mecânico do tratamento inibidor da rapamicina). Em particular, os autores analisaram os perfis de expressão de células TNBC (MDA-MB-231) e epiteliais não malignas da mama (MCF10A) expostas a condições normoxicais (21% O2) e hipóxicas (0,5% O2) e/ou tratadas com um inibidor mTORC1 e -2 PP242. Cada uma das amostras foi sequenciada para mRNA total (T) e polissómico (P). O conjunto de dados contém 32 amostras, representando duas réplicas biológicas para cada combinação de linha celular, nível de oxigênio, status do tratamento e fração do mRNA.
Análise explicativa do conjunto de dados processados em GREIN (Fig. 2) mostra que a mais forte fonte de variação entre as amostras vem das diferenças entre as duas linhas celulares. Isto é reforçado pela análise de correlação dos perfis de expressão completa (Fig. 2A), o agrupamento hierárquico dos 500 genes altamente variáveis com base no desvio absoluto mediano (Fig. 2B), o gráfico PCA 3D das amostras (Fig. 2C) e o gráfico 2D t-SNE (Fig. 2D). Além disso, correlações elevadas entre perfis de expressão para a mesma linha celular (Fig. 2A) indicam bom sinal-ruído nas medições de expressão gênica. A subestrutura adicional de dados indicada pelo gráfico 2D t-SNE foi examinada através de amostras de pintura de acordo com diferentes atributos (Fig. S1 Suplementar). Esta análise revelou que as separações dentro de cada linha celular são induzidas por diferentes frações de mRNA e então diferenças entre condições experimentais.
Next, usamos GREIN para realizar análise estatística de potência com base no padrão de variabilidade biológica observado neste conjunto de dados. Consideramos perfis transcripcionais de cada linha celular exposta à hipóxia e tratada com ou sem PP242, o que leva a quatro comparações. Assumindo uma diferença de expressão de pelo menos duas vezes entre os grupos, na significância estatística de α = 0,01, e com apenas duas réplicas em cada grupo, o poder estatístico de um gene a ser detectado como diferentemente expresso é inferior a 0,55 em todas as comparações (Tabela 2). Nossa análise indica que seriam necessárias quatro réplicas por grupo para se conseguir 80% de potência detectando duas vezes a mudança de expressão (Tabela 2 e Fig. 3A). Em um experimento típico RNA-seq, uma profundidade de seqüenciamento de 20-30 milhões é suficiente para quantificar a expressão gênica para quase todos os genes4,22 o que também é evidente neste conjunto de dados. Também avaliamos o poder estatístico de cada gene a ser detectado como diferentemente expresso no gráfico ‘Detectabilidade dos genes’. O log de contagens médias por milhão (CPM) dos genes foi plotado em relação ao coeficiente de variação biológico dos genes (BCOV) e a potência foi calculada para os genes correspondentes (Fig. 3B). Uma taxa controlada de falsas descobertas de 0,05 e uma percentagem esperada de verdadeiros positivos de 10% foi usada para estimar a significância estatística. Definimos um gene a ser detectável como diferentemente expresso em condição hipóxica se a sua potência for 0,8 ou superior. Como esperado, existe uma relação inversa entre BCOV e potência (Fig. 3B). Além disso, a potência para detectar a expressão diferencial de um gene aumenta com um log CPM ou tamanho de efeito maior.
Um dos objetivos do estudo foi analisar as mudanças transcripcionais em condições hipóxicas e normoxicais com e sem tratamento com PP242 em ambas as linhas de células MCF10A e MDA-MB-231. Criamos assinaturas transcritivas de hipóxia e hipóxia + PP242 no mRNA total por análise de expressão diferencial entre amostras de hipóxia e hipóxia + PP242 respectivamente contra as amostras de controle, enquanto ajustamos para efeito de lote, tratando ‘replicar’ como uma covariada, para cada linha celular separadamente. Encontramos um número maior de genes diferentemente expressos (DE) nas linhas de células MCF10A em comparação ao MDA-MB-231 tanto na hipóxia como na hipóxia + PP242 (Fig. 4A) indicando que talvez a linha de células tumorais esteja melhor equipada para lidar com a hipóxia. Esta análise também mostrou que a maioria dos genes expressos não diferentemente também não são detectáveis, indicando que eles podem representar resultados falsos negativos. Isto está de acordo com a análise de potência mostrando que 4 amostras por grupo seriam necessárias para identificar consistentemente genes expressos diferentemente com BCOV médio. Para identificar genes com menor expressão seria necessário um tamanho de amostra ainda maior.
Expressão diferencial e detectabilidade dos genes. (A) O número de genes (escala log10) não diferentemente expressos e não detectáveis (NDE&NDT), diferentemente expressos (DE), e não diferentemente expressos mas detectáveis (NDE&DT) nas comparações com a normoxia para a fração total de mRNA. Chamamos um gene detectável (DT) se sua potência for ≥0.8 e diferenciável se FDR < 0.05. (B) O gráfico de detectabilidade de genes para a primeira comparação (MCF10A e hipóxia) que visualiza a lista de genes acima mencionada juntamente com suas respectivas alterações de pregas (FC).
Para interpretar genes diferentemente expressos em termos de vias biológicas afetadas, submetemos as assinaturas de expressão diferencial de genes de hipóxia a ferramentas de enriquecimento on-line (DAVID23, ToppGene24, Enrichr25 e Reactome26) via iLINCS. As assinaturas submetidas incluíram uma lista combinada de genes DE e NDE&DT representando prováveis verdadeiros positivos e verdadeiros negativos. Os genes foram selecionados com base em um corte de 0,7 e 0,01 para potência estatística e FDR respectivamente. A Figura 5 ilustra os resultados de enriquecimento obtidos do ToppGene para a assinatura da hipóxia MCF10. Significativamente enriquecida (FDR < 0.05) As categorias top 10 de ontologia genética (GO) do ToppGene e da ferramenta de anotação funcional DAVID incluem resposta à hipóxia, resposta à diminuição dos níveis de oxigênio, angiogênese, regulação da proliferação celular, processo de oxidação-redução e resposta ao estímulo abiótico que são comuns em ambas as linhas celulares (Tabela Complementar S2 e Tabela Complementar S3). A maioria destas categorias é consistente com o estudo original. Além disso, o conjunto ToppGene identificou a rede de fator de transcrição induzido pela hipoxia (HIF-1-alfa) que foi ativada em ambas as linhas celulares (Tabela Suplementar S4 e Tabela Suplementar S5).
Snapshot de algumas das categorias de caminho significativo e ontologia gênica (GO) do ToppGene via iLINCS. Estas categorias são encontradas na comparação entre hipoxia e normoxia na linha celular MCF10A usando uma lista combinada de genes DE e NDE&DT. A linha vertical vermelha é o corte selecionado de 0.05.
Finalmente, nós utilizamos a conexão GREIN com iLINCS para “conectar” a assinatura carregada com assinaturas de derrubada do gene LINCS27 consensus (CGS)18. Encontramos 3.727 assinaturas de derrubada do gene do consenso LINCS que foram significativamente (pValor < 0,05) conectadas com a nossa assinatura carregada. Os genes alvo das 100 principais assinaturas conectadas foram selecionados para uma análise de enriquecimento adicional. Encontramos resposta celular à hipóxia e regulação do Fator induzível pela Hipóxia (HIF) pelo oxigênio na lista dos 10 principais caminhos ativados em ambas as linhas celulares (Tabela Complementar S6 e Tabela Complementar S7). Embora essa análise produza categorias funcionais enriquecidas semelhantes à análise inicial de enriquecimento, ela complementa a análise original, implicando vários genes-alvo que não são expressos de forma diferente, embora sejam suficientemente expressos para serem detectáveis de acordo com nossas análises de potência. A ligação entre estes dois resultados implica estes genes como potenciais reguladores de nível superior da resposta à hipóxia.