Entering edit mode
2.8 years ago
jyotikataria2008
▴
10
I have been working on a 5 umaps which takes annotations from the 5 databases and annotate single cells. My aim is to annotate a cell type with same color across five databases in 5 umaps. I wrote this code -
setwd("F:\NEW\Projects")
pacman::p_load(dplyr, Seurat, patchwork, scSorter, SingleR, celldex,,writexl, ggplot2, stringr, tidyverse)
pbmc.data <- Read10X(data.dir = "filtered_feature_bc_matrix/")
pbmc <- CreateSeuratObject(counts = pbmc.data, project = "pbmc", min.cells = 3, min.features = 200)
pbmc[["percent.mt"]] <- PercentageFeatureSet(pbmc, pattern = "^MT-")
pbmc <- NormalizeData(pbmc, normalization.method = "LogNormalize", scale.factor = 10000)
pbmc <- FindVariableFeatures(pbmc, selection.method = "vst", nfeatures = 2000)
all.genes <- rownames(pbmc)
pbmc <- ScaleData(pbmc, features = all.genes)
pbmc <- RunPCA(pbmc, features = VariableFeatures(object = pbmc))
pbmc <- FindNeighbors(pbmc, dims = 1:10)
pbmc <- FindClusters(pbmc, resolution = 0.5)
Creating different objects for 5 databases
hpca.data <- pbmc
bped.data <- pbmc
dice.data <- pbmc
hd.data <- pbmc
mid.data <- pbmc
Loading data
hpca.se <- HumanPrimaryCellAtlasData()
bped.se <- BlueprintEncodeData()
dice.se <- DatabaseImmuneCellExpressionData()
hd.se <- NovershternHematopoieticData()
mid.se <- MonacoImmuneData()
Setting labels
hpca.se$label.main <- str_replace_all(hpca.se$label.main, c("NK_cell" = "NK cells", "B_cell" = "B cells", "DC" = "Dendritic cells","HSC_-G-CSF" = "HSC G-CSF", "Monocyte" = "Monocytes", "Astrocyte" = "Astrocytes", "Erythroblast" = "Erythroblasts", "Macrophage" = "Macrophages","BM" = "BMs", "MSC"="MSCs","CMP" = "CMPs","GMP" = "GMPs", "MEP"="MEPs", "Myelocyte"="Myelocytes","Neuroepithelial_cell"="Neuroepithelial cells", "T_cells" = "T cells", "iPS_cells" = "iPS cells","Endothelial_cells" = "Endothelial cells", "Tissue_stem_cells" = "Tissue stem cells","Embryonic_stem_cells" = "Embryonic stem cells","Smooth_muscle_cells" = "Smooth muscle cells","Epithelial_cells" = "Epithelial cells"))
hpca.se$label.main[hpca.se$label.main == "HSC_CD34+"] <-"HSC CD34"
hpca.se$label.main[hpca.se$label.main == "Pro-B cells_CD34+"] <-"Pro-B cells CD34 Pos"
hpca.se$label.main[hpca.se$label.main == "Pre-B cells_CD34-"] <-"Pre-B cells CD34 Neg"
bped.se$label.main <- str_replace_all(bped.se$label.main, c("B-cells" = "B cells", "HSC" = "HSCs","DC" = "Dendritic cells"))
bped.se$label.main[bped.se$label.main == "CD4+ T-cells"] <- "CD4 T cells"
bped.se$label.main[bped.se$label.main == "CD8+ T-cells"] <- "CD8 T cells"
dice.se$label.main[dice.se$label.main == "T cells, CD4+"] <- "CD4 T cells"
dice.se$label.main[dice.se$label.main == "T cells, CD8+"] <- "CD8 T cells"
hd.se$label.main[hd.se$label.main == "CD4+ T cells"] <- "CD4 T cells"
hd.se$label.main[hd.se$label.main == "CD8+ T cells"] <- "CD8 T cells"
mid.se$label.main[mid.se$label.main == "CD4+ T cells"] <- "CD4 T cells"
mid.se$label.main[mid.se$label.main == "CD8+ T cells"] <- "CD8 T cells"
#------Humanprimarycellatlas
pred.hesc <- SingleR(test = counts, ref = hpca.se, assay.type.test=1,
labels = hpca.se$label.main)
pred.hesc.transform <- t(pred.hesc)
hesc_cells_names <- pred.hesc.transform$labels
hpca.data <- AddMetaData(object = hpca.data,metadata = hesc_cells_names,col.name = 'type')
Idents(hpca.data) <- hpca.data@meta.data$type
##---- BlueprintEncodeData()
pred.bped <- SingleR(test = counts, ref = bped.se, assay.type.test=1,
labels = bped.se$label.main)
pred.bped.transform <- t(pred.bped)
bped_cells_names <- pred.bped.transform$labels
bped.data <- AddMetaData(object = bped.data,metadata = bped_cells_names,col.name = 'type')
Idents(bped.data) <- bped.data@meta.data$type
DatabaseImmuneCellExpressionData()
pred.dice <- SingleR(test = counts, ref = dice.se, assay.type.test=1,
labels = dice.se$label.main)
pred.dice.transform <- t(pred.dice)
dice_cells_names <- pred.dice.transform$labels
dice.data <- AddMetaData(object = dice.data,metadata = dice_cells_names,col.name = 'type')
Idents(dice.data) <- dice.data@meta.data$type
NovershternHematopoieticData()
pred.hd <- SingleR(test = counts, ref = hd.se, assay.type.test=1,
labels = hd.se$label.main)
pred.hd.transform <- t(pred.hd)
hd_cells_names <- pred.hd.transform$labels
hd.data <- AddMetaData(object = hd.data,metadata = hd_cells_names,col.name = 'type')
Idents(hd.data) <- hd.data@meta.data$type
MonacoImmuneData()
pred.mid <- SingleR(test = counts, ref = mid.se, assay.type.test=1,
labels = mid.se$label.main)
pred.mid.transform <- t(pred.mid)
mid_cells_names <- pred.mid.transform$labels
mid.data <- AddMetaData(object = mid.data,metadata = mid_cells_names,col.name = 'type')
Idents(mid.data) <- mid.data@meta.data$type
colors50<- c("#a6cee3","#1f78b4","#b2df8a","#33a02c","#fb9a99","#e31a1c","#fdbf6f","#ff7f00","#cab2d6","#6a3d9a","#ffff99","#b15928","#01665e","#c51b7d","#67001f","#53061","#4d4d4d","#f768a1","#276419","#8e0152","#7a0177","#80cdc1","#f1b6da","#e0e0e0","#f46d43","#878787","#045a8d","#fcc5c0","#1a1a1a","#4d004b","#4d4d4d","#f768a1","#276419","#8e0152","#7a0177","#80cdc1","#f1b6da","#e0e0e0","#f46d43","#878787","#045a8d","#fcc5c0","#1a1a1a","#4d004b","#4d4d4d","#f768a1","#276419","#8e0152","#7a0177","#80cdc1","#f1b6da","#e0e0e0","#f46d43","#878787","#045a8d","#fcc5c0","#1a1a1a","#4d004b")
hpca.data <- RunUMAP(hpca.data, dims = 1:10)
bped.data <- RunUMAP(bped.data, dims = 1:10)
dice.data <- RunUMAP(dice.data, dims = 1:10)
hd.data <- RunUMAP(hd.data, dims = 1:10)
mid.data <- RunUMAP(mid.data, dims = 1:10)
hpca.clusters <- Idents(hpca.data)
bped.clusters <- Idents(bped.data)
dice.clusters <- Idents(dice.data)
hd.clusters <- Idents(hd.data)
mid.clusters <- Idents(mid.data)
Assigning same colors to the same cell annotated across five databases
clusters <- c(hpca.clusters,bped.clusters,dice.clusters,hd.clusters,mid.clusters)
unique.clusters <- unique(clusters)
cluster.colors <- colors50[1:length(unique.clusters)]
names(cluster.colors) <- unique.clusters
cols <- cluster.colors[order(as.integer(names(cluster.colors)))]
DimPlot(object = hpca.data, reduction = "umap", cols = cluster.colors)
DimPlot(object = hpca.data, reduction = "umap", cols = cols)
These both gave me grey umap in my rmd on the server but its showing me correct colours with the rstudio. Is the code correct and well logically made?
