Markers identified by Seurat:
## see code_processing.md
## FindAllMarkers was run on the integrated Seurat object
srt_integrated <- FindClusters(srt_integrated, resolution = .3)
markers_seurClsts.3 <- FindAllMarkers(srt_integrated)
## add cluster information back to the SingleCellExperiment object
cd_srt <- srt_integrated@meta.data
scf$integClust_res.3 <- cd_srt[colnames(scf),]$integrated_snn_res.0.3
Markers identified with scran:
label_de_gns_ctrl <- scran::findMarkers(scf[, scf$condition == "Ctrl"],
as.character(scf[,scf$condition == "Ctrl"]$labels),
assay.type = "log1p_sctransform",
direction = "up", lfc = 0)
## extracting top hits
ldg <- scABC2::extract_markers(scf, label_de_gns, FDR_thresh = 0.01, rank_thresh = 10)
ldg.un <- ldg[!grepl("^RP[LS]", gene_symbol) & classify == "unique"]
The development version of SingleR was used.
library(SingleR) # v.0.99
hpca.se <- HumanPrimaryCellAtlasData()
common <- intersect(rownames(scf), rownames(hpca.se))
hpca.se <- hpca.se[common,]
scf1 <- scf[common,]
pred_hpca <- SingleR(test = scf1, ref = hpca.se, labels = hpca.se$label.main)
save(pred_hpca, file = "SingleR_labels_HPCA.rda")
Downloaded matrix of TPM (not log2-transformed) from here. The labels come from the supplemental tables of Nowakowski et al., 2017. Spatiotemporal gene expression trajectories reveal developmental hierarchies of the human cortex..
## matrix of TPM (not log-transformed) ------------------------------------------
## https://cells.ucsc.edu/cortex-dev/exprMatrix.tsv.gz
exprs.nowa <- fread(cmd=paste0("gunzip -c ", dfold, "exprMatrix.tsv.gz"))
logtpm <- log1p(as.matrix(exprs.nowa[, -"gene", with=FALSE]))
rownames(logtpm) <- exprs.nowa$gene
### colData: cell labels -------------------------------------------------------
cellinf <- openxlsx::read.xlsx(paste0(dfold, "aap8809_Nowakowski_SM-Tables-S1-S11.xlsx"), sheet = 3) %>%
as.data.table
celltypeinf <- openxlsx::read.xlsx(paste0(dfold,"aap8809_Nowakowski_SM-Tables-S1-S11.xlsx"), sheet = 4) %>%
as.data.table
setnames(celltypeinf, "Cluster.Name", "WGCNAcluster")
cellinf <- celltypeinf[, -2, with=FALSE] %>% .[cellinf, on = "WGCNAcluster"]
cellinf <- cellinf[ !grepl("^Unknown", Cluster.Interpretation)]
cold <- DataFrame(as.data.frame(cellinf))
rownames(cold) <- cold$Cell
## create SE -------------------------------------------------------------------
nowa.sce <- SingleCellExperiment(assays = list(logcounts = logtpm),
colData = cold[colnames(logtpm),])
nowa.sce <- nowa.sce[, !is.na(nowa.sce$WGCNAcluster)]
common <- intersect(rownames(nowa.sce), rownames(scf))
scf2 <- scf[common,]
nowa.sce <- nowa.sce[common,]
## aggregate the logTPMs --------------------------------------------------
nowa.aggr <- aggregateReference(nowa.sce,
labels = nowa.sce$WGCNAcluster,
power=0.5)
## run SingleR------------------------------------------------------------
pred_nowa_aggregated <- SingleR(test = scf2, ref = nowa.aggr,
labels = nowa.aggr$label,
assay.type.test = "log1p_sctransform")
Based on the insights from the marker gene analyses as well as the automated cell annotation, we labelled individual cells as follows:
- cluster 4: neurons
- clusters 1, 2, 7: progenitors
- clusters 3, 6: proliferating
- cluster 9: endothelial
- cluster 10: astrocytes
- clusters 0, 11, 12,: radial glia (RG)
- cluster 8: myeloid
- cluster 5: mural
scf$labels <- ifelse(scf$integClust_res.3 == "4", "Neurons",
ifelse(scf$integClust_res.3 %in% c("2","7","1"), "Progs.",
ifelse(scf$integClust_res.3 %in% c("3","6"), "Proliferating",
ifelse(scf$integClust_res.3 == "9", "Endothelial",
ifelse(scf$integClust_res.3 == "10", "Astrocytes",
ifelse(scf$integClust_res.3 %in% c("0","12","11"), "RG",
ifelse(scf$integClust_res.3 == "8", "Myeloid",
ifelse(scf$integClust_res.3 == "5", "Mural", NA))))))))