1 Visão geral

Esta vinheta demonstra como estimar misturas de tipos celulares na seção anterior do cérebro de camundongo stxBrain com o método Robust Cell Type Decomposition (RCTD), usando o Seurat 5.5.1 e a versão 2.2.1 do pacote spacexr disponível no GitHub (Cable et al. 2022).

O objeto espacial será obtido diretamente pelo pacote SeuratData, usando a seção anterior1 (Hao et al. 2024). Nenhum objeto espacial previamente processado será reutilizado. Essa separação é importante porque o RCTD modela contagens brutas de expressão e, por padrão, exige valores inteiros. Matrizes normalizadas, como as produzidas por NormalizeData() ou SCTransform(), não devem ser fornecidas como entrada do método.

Ao final do roteiro, você será capaz de:

  1. explicar por que a deconvolução requer uma referência single-cell anotada;
  2. obter do SeuratData as contagens brutas e as coordenadas do Visium;
  3. preparar a referência single-cell do córtex de camundongo do Allen Institute (Tasic et al. 2016);
  4. executar o RCTD;
  5. inspecionar e visualizar as proporções estimadas; e
  6. acrescentar ao objeto Seurat a classificação e os tipos celulares atribuídos pelo RCTD.

2 Fundamentos conceituais

Um spot do Visium é maior que uma única célula e pode capturar RNA proveniente de vários tipos celulares. O RCTD aprende perfis de expressão a partir de dados anotados de scRNA-seq (single-cell) ou snRNA-seq (single-nucleus), corrige diferenças entre as plataformas da referência e da amostra espacial e estima a contribuição de cada tipo celular em cada spot (Cable et al. 2022).

A análise possui duas entradas e uma saída principal:

Componente Objeto utilizado Informação necessária
Amostra espacial stxBrain, seção anterior1, carregada por SeuratData contagens brutas de genes por spots e coordenadas dos spots
Referência single-cell córtex de camundongo do Allen Institute contagens brutas de genes por células e um rótulo de tipo celular por célula
Resultado do RCTD objeto da classe RCTD classificação e pesos de até dois tipos celulares por spot

O spacexr oferece três modos:

  • doublet: ajusta no máximo dois tipos celulares por spot;
  • multi: ajusta até MAX_MULTI_TYPES tipos celulares por spot; e
  • full: não impõe um limite fixo ao número de tipos celulares.

Usaremos o modo doublet, por ser mais leve. No modo doublet, o RCTD atribui um ou dois tipos celulares a cada spot. Isso é uma simplificação do modelo e não significa que um spot Visium contenha, biologicamente, no máximo duas células ou dois tipos celulares.

3 Antes de começar

3.1 Requisitos computacionais

  • R e RStudio;
  • no mínimo 20 GB de memória RAM, embora uma quantidade maior seja recomendada;
  • acesso à internet para a primeira instalação dos pacotes e dos dados; e
  • múltiplos núcleos de CPU para reduzir o tempo de execução, quando disponíveis.

Como referência de desempenho, a execução validada para a etapa de cálculo do algortimo do RCTD levou 14,22 minutos em um servidor de alta performance com 160 GB de RAM e 30 núcleos de CPU disponíveis para processamento paralelo. O tempo pode variar bastante de acordo com o equipamento, a carga do servidor e o número de tipos celulares da referência. O parâmetro max_cores controla o paralelismo na CPU; esta análise não utiliza aceleração por GPU.

Além disso, a paralelização empregada neste fluxo está disponível apenas em Linux e macOS. No Windows, max_cores deve ser igual a 1, o que pode aumentar consideravelmente o tempo de processamento. Diferentemente das demais, esta vinheta inclui figuras previamente geradas para que o fluxo possa ser acompanhado mesmo quando não for possível executar a análise completa durante o curso.

Durante o curso, recomenda-se executar os blocos de código na ordem em que aparecem. Assim, cada objeto intermediário pode ser conferido antes de prosseguir.

3.2 Estrutura esperada das pastas

O conjunto stxBrain não precisa ser copiado para o repositório: ele é instalado pelo SeuratData na biblioteca de pacotes do R. O único arquivo de entrada local é allen_cortex.rds, que será baixado mais adiante e armazenado em data/. Os resultados serão gravados em resultados/ e as figuras em imagens/.

Curso_ST/
|-- data/
|   `-- allen_cortex.rds
|-- imagens/
`-- resultados/

4 Conferir os pacotes necessários

A instalação de todos os pacotes e de suas dependências está centralizada em Configuracao_do_ambiente_RStudio.Rmd. O roteiro pressupõe esse ambiente recém-instalado e adota o Seurat 5.5.1 como referência, sem caminhos alternativos para versões anteriores. Para a deconvolução, utiliza a versão 2.2.1 do repositório original do spacexr no GitHub, cuja interface oferece create.RCTD() e run.RCTD() e não é intercambiável com a interface do Bioconductor.

pacotes_rctd <- c(
  "Seurat", "SeuratObject", "SeuratData", "Matrix", "ggplot2", "spacexr"
)

pacotes_rctd_ausentes <- pacotes_rctd[
  !vapply(pacotes_rctd, requireNamespace, logical(1), quietly = TRUE)
]

if (length(pacotes_rctd_ausentes) > 0) {
  stop(
    "Execute primeiro Configuracao_do_ambiente_RStudio.Rmd. ",
    "Pacotes ausentes: ",
    paste(pacotes_rctd_ausentes, collapse = ", ")
  )
}

5 Definir os diretórios de trabalho

Para que cada bloco de código possa ser executado separadamente, o fluxo usa a pasta de trabalho atual como ponto de partida. Antes de continuar, confirme que ela corresponde à pasta vinhetas/ do curso.

Siga estes passos:

  1. Confirme que o diretório de trabalho é a pasta que contém esse arquivo. Se necessário, use o menu Session > Set Working Directory > To Source File Location.
  2. Execute getwd() no Console. O R deverá mostrar um caminho terminado em Curso_ST/vinhetas.
  3. Caso o caminho seja diferente, ajuste o diretório de trabalho antes de executar o bloco abaixo.
diretorio_vinhetas <- getwd()

diretorio_curso <- dirname(diretorio_vinhetas)
diretorio_dados <- file.path(diretorio_curso, "data")
diretorio_imagens <- file.path(diretorio_curso, "imagens")
diretorio_resultados <- file.path(diretorio_curso, "resultados")

dir.create(diretorio_dados, showWarnings = FALSE, recursive = TRUE)
dir.create(diretorio_imagens, showWarnings = FALSE, recursive = TRUE)
dir.create(diretorio_resultados, showWarnings = FALSE, recursive = TRUE)

diretorio_curso
diretorio_dados
diretorio_imagens
diretorio_resultados

6 Carregar os pacotes

library("Seurat")
library("SeuratData")
library("Matrix")
library("ggplot2")
library("spacexr")

options(stringsAsFactors = FALSE)

Confira as versões carregadas. Os resultados esperados são Seurat 5.5.1 e spacexr 2.2.1 ou maior:

packageVersion("Seurat")
packageVersion("spacexr")

7 Obter o stxBrain diretamente do SeuratData

7.1 Conferir o conjunto de dados

Nesta vinheta, o objeto Seurat stxBrain será carregado novamente porque o RCTD precisa das contagens brutas, e não do objeto normalizado produzido nas etapas anteriores do curso. O bloco abaixo apenas confirma que o conjunto de dados já está instalado. Caso ele ainda não esteja disponível, retorne à vinheta Configuracao_do_ambiente_RStudio.Rmd.

if (!requireNamespace("stxBrain.SeuratData", quietly = TRUE)) {
  stop(
    "O conjunto stxBrain não está instalado.\nExecute primeiro a vinheta ",
    "Configuracao_do_ambiente_RStudio.Rmd: \n",
    "https://flaviogckessler.github.io/Curso_ST/vinhetas/Configuracao_do_ambiente_RStudio.html"
  )
}

message("O conjunto stxBrain já foi instalado.")

7.2 Carregar a seção anterior 1

Com o objeto instalado, basta realizar o seu carregamento:

stxbrain <- SeuratData::LoadData("stxBrain", type = "anterior1")

stxbrain

7.3 Extrair as contagens espaciais brutas

O RCTD receberá somente a camada counts do ensaio Spatial, que conserva as contagens inteiras e não normalizadas. As camadas data e scale.data, assim como um ensaio SCT, não devem ser usadas como entrada.

contagens_espaciais <- SeuratObject::LayerData(
  object = stxbrain,
  assay = "Spatial",
  layer = "counts"
)

print(paste("Número de spots:",dim(contagens_espaciais)[2]))
print(paste("Número de genes:",dim(contagens_espaciais)[1]))

contagens_espaciais[
   seq_len(min(5L, nrow(contagens_espaciais))),
   seq_len(min(5L, ncol(contagens_espaciais))),
   drop = FALSE]

7.4 Examinar a profundidade de sequenciamento

O total de transcritos varia entre os spots. Este resumo permite examinar a profundidade de sequenciamento e a ordem de grandeza das bibliotecas antes da deconvolução.

summary(Matrix::colSums(contagens_espaciais))

7.5 Extrair as coordenadas dos spots

No Seurat 5.5.1, a imagem espacial do objeto stxBrain fornece as coordenadas não escaladas nas colunas x e y. O identificador de cada spot permanece nos nomes das linhas, permitindo alinhar diretamente a tabela de coordenadas às colunas da matriz de contagens antes de criar o objeto do spacexr.

nome_imagem <- SeuratObject::Images(stxbrain)[1]
coordenadas_espaciais <- SeuratObject::GetTissueCoordinates(
  object = stxbrain,
  image = nome_imagem,
  scale = NULL
)

coordenadas_espaciais <- as.data.frame(
  coordenadas_espaciais[, c("x", "y"), drop = FALSE]
)
coordenadas_espaciais <- coordenadas_espaciais[
  colnames(contagens_espaciais), , drop = FALSE
]

head(coordenadas_espaciais)

8 Preparar a referência do córtex do Allen Institute

O RCTD requer uma referência de RNA-seq biologicamente adequada e anotada. Usaremos a taxonomia celular do córtex adulto de camundongo do Allen Institute, também empregada no tutorial espacial do Seurat (Tasic et al. 2016). Como a referência e a amostra espacial são de camundongo, os genes podem ser associados por símbolo gênico.

8.1 Baixar a referência uma única vez

Baixe a referência single-cell com o código abaixo. O arquivo será salvo como data/allen_cortex.rds. Atenção: o download possui cerca de 1 GB e precisa ser realizado apenas uma vez.

arquivo_referencia <- file.path(diretorio_dados, "allen_cortex.rds")

download.file(
  url = paste0(
    "https://www.dropbox.com/scl/fi/sgmfzk72cdoc945ap19c9/",
    "allen_cortex.rds?rlkey=qfou47stxczo14m37vmea54wl&dl=1"
  ),
  destfile = arquivo_referencia,
  mode = "wb"
)

8.2 Extrair as contagens e os rótulos da referência

A referência publicada armazena os rótulos de tipos celulares na coluna subclass.

Cuidado: o carregamento do arquivo de referência, células corticais de camundongo disponibilizadas pelo Allen Institute, é realizado no código abaixo. O arquivo é grande e vai ocupar 3GB de memória quando for carregado para o ambiente de trabalho.

arquivo_referencia <- file.path(diretorio_dados, "allen_cortex.rds")

# Confere se arquivo existe
if (!file.exists(arquivo_referencia)) {
  stop(
    "A referência single-cell não foi encontrada em: ",
    arquivo_referencia,
    "\nBaixe o arquivo como descrito acima ou copie o arquivo fornecido pelo instrutor para data/."
  )
}

# Carrega o arquivo de referência
referencia_allen <- readRDS(arquivo_referencia)

#Extrai metadados das células do arquivo de referência
metadados_referencia <- referencia_allen[[]]

# Extrai contagens do arquivo de referência
contagens_referencia <- SeuratObject::LayerData(
  object = referencia_allen,
  assay = "RNA",
  layer = "counts"
)
head(contagens_referencia[,1:5])

# Organiza os metadados para seguir a mesma ordem da matriz de contagem
metadados_referencia <- metadados_referencia[
  colnames(contagens_referencia), , drop = FALSE
]

head(metadados_referencia)
# Extrai coluna dos tipos celulares dos metadados
tipos_celulares_referencia <- as.character(
  metadados_referencia$subclass
)
# Corrige strings dos nomes das células
tipos_celulares_referencia <- gsub(
  "/", "-", tipos_celulares_referencia, fixed = TRUE
)

celulas_validas <- !is.na(tipos_celulares_referencia) &
  nzchar(tipos_celulares_referencia)

tamanhos_tipos <- table(tipos_celulares_referencia[celulas_validas])

# Filtra tipos celulares que estão presentes em menos de 25 células
tipos_elegiveis <- names(tamanhos_tipos[tamanhos_tipos >= 25])
celulas_validas <- celulas_validas &
  tipos_celulares_referencia %in% tipos_elegiveis

# Retira células com tipo celular raro da contagem de referência
contagens_referencia <- contagens_referencia[
  , celulas_validas, drop = FALSE
]
tipos_celulares_referencia <- droplevels(
  factor(tipos_celulares_referencia[celulas_validas])
)

print(paste("N° de células de referência: ",ncol(contagens_referencia)))
print("Células de referência:")
levels(tipos_celulares_referencia)

O filtro padrão do RCTD exige pelo menos 25 células por tipo celular. Os tipos raros são removidos explicitamente acima, evitando que desapareçam de forma silenciosa durante o pré-processamento.

8.3 Identificar genes compartilhados

Vamos analisar os dados gerais que iremos trabalhar

genes_comuns <- intersect(
  rownames(contagens_espaciais),
  rownames(contagens_referencia)
)

resumo_genes <- data.frame(
  fonte = c("amostra espacial", "referência single-cell", "interseção"),
  numero_de_genes = c(
    nrow(contagens_espaciais),
    nrow(contagens_referencia),
    length(genes_comuns)
  )
)
resumo_genes

if (length(genes_comuns) < 1000L) {
  stop(
    "Foram encontrados menos de 1.000 genes compartilhados. Verifique o organismo, ",
    "os identificadores gênicos e se ambas as matrizes contêm contagens brutas."
  )
}

contagens_espaciais <- contagens_espaciais[
  genes_comuns, , drop = FALSE
]
contagens_referencia <- contagens_referencia[
  genes_comuns, , drop = FALSE
]

# Retirar dados pesados do ambiente de trabalho (liberar RAM)
rm(resumo_genes, referencia_allen)
invisible(gc())

O dado espacial (stxBrain) possui 31053 genes identificados. A referência single-cell com os tipos celulares anotados possui 34617 genes identificados. O dado espacial e a referência possuem uma intersecção de 19293 genes identificados.

9 Construir os objetos de entrada do spacexr

O spacexr recebe a referência single-cell em um objeto Reference e a amostra espacial em um objeto SpatialRNA. Em ambos, as somas por coluna representam o número total de UMIs de cada célula ou spot.

names(tipos_celulares_referencia) <- colnames(contagens_referencia)

# Construir objeto do pacote spacexr para o dado da referência
referencia_rctd <- spacexr::Reference(
  contagens_referencia,
  tipos_celulares_referencia,
  Matrix::colSums(contagens_referencia)
)

# Construir objeto do pacote spacexr para o dado espacial
consulta_rctd <- spacexr::SpatialRNA(
  as.data.frame(coordenadas_espaciais),
  contagens_espaciais,
  Matrix::colSums(contagens_espaciais)
)

data.frame(objeto = c("referencia_rctd", "consulta_rctd"),
  genes = c(nrow(referencia_rctd@counts),nrow(consulta_rctd@counts)),
  unidades = c(ncol(referencia_rctd@counts),ncol(consulta_rctd@counts)),
  tipo_de_unidade = c("células", "spots"))

# Retirar dados pesados do ambiente de trabalho (liberar RAM)
rm(contagens_espaciais,contagens_referencia,coordenadas_espaciais,
  metadados_referencia,tipos_celulares_referencia,genes_comuns)
invisible(gc())

10 Executar o RCTD

10.1 Criar o objeto RCTD

create.RCTD() prepara os dados, seleciona genes informativos e aprende os perfis dos tipos celulares da referência. No código abaixo foi atribuído 1 núcleo de processamento para evitar exaurir o seu computador. Contudo, isso irá deixar a tarefa muito demorado. Consulte o número de cores que o CPU do seu computador possui e assine um número que seja seguro. A parelelização está apenas disponíveis em sistemas Linux e macOS. No Windows, use obrigatoriamente 1.

set.seed(20260904)

numero_nucleos <- 1

RCTD <- spacexr::create.RCTD(
  consulta_rctd,
  referencia_rctd,
  max_cores = numero_nucleos
)

rm(consulta_rctd, referencia_rctd)
invisible(gc())

10.2 Executar a deconvolução

O modo doublet procura um ou dois tipos celulares por spot. O tempo de execução varia de acordo com o servidor e o número de células da referência.

Cuidado: Essa é a etapa mais demorada da análise!

inicio_rctd <- Sys.time()

RCTD <- spacexr::run.RCTD(
  RCTD,
  doublet_mode = "doublet"
)

fim_rctd <- Sys.time()
tempo_rctd <- fim_rctd - inicio_rctd
tempo_rctd

saveRDS(
  RCTD,
  file = file.path(
    diretorio_resultados,
    "stxbrain_RCTD_doublet_results.rds"
  )
)

rm(inicio_rctd, fim_rctd, tempo_rctd, numero_nucleos)
invisible(gc())

11 Inspecionar o resultado do RCTD

No modo doublet, os resultados principais ficam em dois componentes do objeto RCTD:

  • results_df: classificação do spot e tipos celulares atribuídos;
  • weights_doublet: pesos associados ao primeiro e ao segundo tipo.
resultado_spots <- RCTD@results$results_df
pesos_doublet <- as.data.frame(RCTD@results$weights_doublet)
colnames(pesos_doublet) <- c("peso_primeiro_tipo", "peso_segundo_tipo")

resumo_spots_rctd <- cbind(resultado_spots, pesos_doublet)

head(resumo_spots_rctd)
table(resumo_spots_rctd$spot_class)

As classes produzidas pelo RCTD são:

  • singlet: um tipo celular foi atribuído com confiança;
  • doublet_certain: dois tipos celulares foram atribuídos com confiança;
  • doublet_uncertain: somente o primeiro tipo é confiável; e
  • reject: o método não produziu uma atribuição confiável.

12 Acrescentar as anotações ao objeto Seurat

Vamos acrescentar a tabela do RCTD ao metadado do stxbrain. O Seurat utiliza os códigos de barras presentes nos nomes das linhas para posicionar cada resultado no spot correto.

stxbrain <- SeuratObject::AddMetaData(
  object = stxbrain,
  metadata = resumo_spots_rctd
)

rm(RCTD, resultado_spots, pesos_doublet)
invisible(gc())

13 Visualizar as classificações no tecido

Primeiro, observe onde o RCTD identificou singlet, doublet_certain, doublet_uncertain ou reject. Ao executar os blocos manualmente, todas as figuras serão salvas como PNG em imagens/, com resolução moderada de 120 dpi.

p_classes <- Seurat::SpatialDimPlot(
  object = stxbrain,
  images = nome_imagem,
  group.by = "spot_class",
  pt.size.factor = 1.6,
  crop = TRUE
) +
  labs(title = "Classificação dos spots pelo RCTD")

p_classes

ggplot2::ggsave(
  filename = file.path(
    diretorio_imagens,
    "rctd_classificacao_spots.png"
  ),
  plot = p_classes,
  width = 7,
  height = 5,
  units = "in",
  dpi = 120,
  bg = "white"
)

rm(p_classes)
invisible(gc())

A figura abaixo foi produzida durante a execução manual validada da vinheta e é apresentada aqui como exemplo. Ela é incorporada ao HTML a partir do arquivo PNG já salvo, sem executar novamente o código.

Classificação dos spots pelo RCTD.
Classificação dos spots pelo RCTD.

Em seguida, visualize o primeiro e o segundo tipo celular atribuídos. O segundo tipo deve ser interpretado somente junto com spot_class: ele é confiável quando o spot foi classificado como doublet_certain.

p1 <- Seurat::SpatialDimPlot(
  object = stxbrain,
  images = nome_imagem,
  group.by = "first_type",
  pt.size.factor = 1.6,
  crop = TRUE
) +
  labs(title = "Primeiro tipo celular")

p2 <- Seurat::SpatialDimPlot(
  object = stxbrain,
  images = nome_imagem,
  group.by = "second_type",
  pt.size.factor = 1.6,
  crop = TRUE
) +
  labs(title = "Segundo tipo celular")

p_tipos <- p1 + p2
p_tipos

ggplot2::ggsave(
  filename = file.path(
    diretorio_imagens,
    "rctd_primeiro_tipo_celular.png"
  ),
  plot = p1,
  width = 7,
  height = 5,
  units = "in",
  dpi = 120,
  bg = "white"
)

ggplot2::ggsave(
  filename = file.path(
    diretorio_imagens,
    "rctd_segundo_tipo_celular.png"
  ),
  plot = p2,
  width = 7,
  height = 5,
  units = "in",
  dpi = 120,
  bg = "white"
)

ggplot2::ggsave(
  filename = file.path(
    diretorio_imagens,
    "rctd_tipos_celulares_combinados.png"
  ),
  plot = p_tipos,
  width = 10,
  height = 5,
  units = "in",
  dpi = 120,
  bg = "white"
)

rm(p1, p2, p_tipos)
invisible(gc())

Os mapas individuais permitem examinar com mais detalhe a distribuição espacial de cada atribuição.

Primeiro tipo celular atribuído a cada spot.
Primeiro tipo celular atribuído a cada spot.
Segundo tipo celular atribuído a cada spot.
Segundo tipo celular atribuído a cada spot.

O painel combinado facilita a comparação direta entre o primeiro e o segundo tipo celular.

Comparação entre o primeiro e o segundo tipo celular atribuídos pelo RCTD.
Comparação entre o primeiro e o segundo tipo celular atribuídos pelo RCTD.

14 Salvar os resultados reutilizáveis

write.csv(
  resumo_spots_rctd,
  file = file.path(
    diretorio_resultados,
    "stxbrain_RCTD_doublet_spot_summary.csv"
  ),
  row.names = TRUE
)

saveRDS(
  stxbrain,
  file = file.path(
    diretorio_resultados,
    "stxbrain_Seurat_com_RCTD.rds"
  )
)

Os arquivos possuem finalidades complementares:

  • stxbrain_RCTD_doublet_results.rds: objeto completo produzido pelo RCTD;
  • stxbrain_RCTD_doublet_spot_summary.csv: classificação, primeiro e segundo tipos celulares e seus pesos por spot; e
  • stxbrain_Seurat_com_RCTD.rds: objeto Seurat com os resultados do RCTD no metadado.

As figuras rctd_classificacao_spots.png, rctd_primeiro_tipo_celular.png, rctd_segundo_tipo_celular.png e rctd_tipos_celulares_combinados.png serão gravadas em imagens/.

15 Lista de verificação para interpretação

Antes de formular conclusões biológicas, verifique:

  1. Qual foi a classificação do spot? Não interprete spots classificados como reject como atribuições confiáveis.
  2. O segundo tipo é confiável? Em doublet_uncertain, somente o primeiro tipo deve ser tratado como confiável.
  3. A referência é adequada? A referência do Allen Institute enfatiza populações corticais; tipos ausentes da referência não podem ser descobertos pelo RCTD.
  4. As previsões são espacialmente plausíveis? Populações esperadas no córtex e na substância branca devem formar padrões anatômicos coerentes, em vez de pontos dispersos aleatoriamente.

As atribuições e os pesos do RCTD são estimativas baseadas em um modelo. Eles não correspondem diretamente ao número de células presentes no spot.

16 Informações da sessão

sessionInfo()

Referências

Cable, Dylan M., Evan Murray, Luli S. Zou, Aleksandr Goeva, Evan Z. Macosko, Fei Chen, and Rafael A. Irizarry. 2022. “Robust Decomposition of Cell Type Mixtures in Spatial Transcriptomics.” Nature Biotechnology 40 (4): 517–26. https://doi.org/10.1038/s41587-021-00830-w.
Hao, Yuhan, Tim Stuart, Madeline H. Kowalski, Saket Choudhary, Paul Hoffman, Austin Hartman, Avi Srivastava, et al. 2024. “Dictionary Learning for Integrative, Multimodal and Scalable Single-Cell Analysis.” Nature Biotechnology 42: 293–304. https://doi.org/10.1038/s41587-023-01767-y.
Tasic, Bosiljka, Vilas Menon, Thuc Nghi Nguyen, Tae Kyung Kim, Tim Jarsky, Zizhen Yao, Boaz Levi, et al. 2016. “Adult Mouse Cortical Cell Taxonomy Revealed by Single Cell Transcriptomics.” Nature Neuroscience 19: 335–46. https://doi.org/10.1038/nn.4216.