Setup

Notes

  • Cell suspensions of meningeal dura enriched by MACS for cells expressing Cd45, Cd31, and Lyve-1
    • Cd45/Ptprc/B220/Ly-5/Lyt-4/T200: transmembrane protein tyrosine phosphatase, located on most haematopoietic cells, positive selection of leukocytes
    • Cd31/Pecam1: adhesion and signaling receptor that is expressed on endothelial and hematopoietic cells, positive selection of endothelial cells
    • Lyve1/1200012G08Rik/Lyve-1/Xlkd1/lymphatic vessel endothelial HA receptor-1: found primarily on lymphatic endothelial cells
  • The meningeal layer of the dura mater is a durable, dense fibrous membrane that passes through the foramen magnum and is continuous with the dura mater of the spinal cord. The meningeal layer of the dura mater creates several dural folds that divide the cranial cavity into freely communicating spaces.
  • Aimed for ~5000 cells per sample to be sequenced at least 50,000 reads per cell
    • Might ended up with more because we sequenced with the NovaSeq S2
  • Expecting to see a lot of leukocytes and vascular cells
  • The purpose is to check the effect of sex, Apoe isoform expression and age (and respective interactions) on the transcriptomes of the different cell subpopulations
  • 8 experimental groups/samples
    • E3_2M_F = APOE3 expressing; 2 months old; female
    • E3_14M_F = APOE3 expressing; 14 months old; female
    • E3_2M_M = APOE3 expressing; 2 months old; male
    • E3_14M_M = APOE3 expressing; 14 months old; male
    • E4_2M_F = APOE4 expressing; 2 months old; female
    • E4_14M_F = APOE4 expressing; 14 months old; female
    • E4_2M_M = APOE4 expressing; 2 months old; male
    • E4_14M_M = APOE4 expressing; 14 months old; male

Working directory

knitr::opts_knit$set(root.dir = ".")

Libraries

# load packages
library(dplyr)        # ungroup()
library(ggtree)       # BuildClusterTree()
library(gridExtra)    # grid.arrange()
library(gtools)       # smartbind()
library(parallel)     # detectCores()
library(plotly)       # plot_ly()
library(Seurat)       # Idents()
library(ShinyCell)    # createConfig()
library(tidyr)        # %>%

# work in parallel
options(mc.cores = detectCores() - 1)

Variables and functions

# variables
sample_order <- c("E3.M","E3.F","E4.M","E4.F")
sample_colors <- c("#26946A","#1814A1","#EAC941","#DF5F00")
sample_order2 <- c("Male E3", "Male E4", "Female E3", "Female E4")
isoform_order <- c("E4","E3")
isoform_colors <- c("darkgray","cornflowerblue")
sex_order <- c("Male","Female")
sex_colors <- c("darkgray","purple")

# single cell functions
source("../../refs/Kennedi_single_cell_functions.R")

# save function
saveToPDF <- function(...) {
    d = dev.copy(pdf,...)
    dev.off(d)
}

Read data

mouse.unannotated <- readRDS("../../rObjects/pass1_unannotated.rds")
Idents(mouse.unannotated) <- "seurat_clusters"
DefaultAssay(mouse.unannotated) <- "RNA"
mouse.unannotated <- NormalizeData(mouse.unannotated)

mouse.unannotated
## An object of class Seurat 
## 42466 features across 21237 samples within 2 assays 
## Active assay: RNA (21389 features, 0 variable features)
##  1 other assay present: SCT
##  2 dimensional reductions calculated: pca, umap

Unannotated QC

UMAP

cluster_colors <- c("black","gray40","gray","red1","blue","magenta1","darkorange2",
                    "darkorange4","yellow1","yellow4","yellow2","green","lightgreen",
                    "chartreuse1","Aquamarine","cyan","SteelBlue","red4","forestgreen",
                    "purple1","purple4","orange","plum1","salmon","tan","chocolate4")
u1 <- DimPlot(object = mouse.unannotated,
              reduction = "umap",
              shuffle = TRUE,
              repel = TRUE,
              cols = cluster_colors,
              label = TRUE)
u1

Feature plots

# UMAP percent.mt
FeaturePlot(mouse.unannotated,
            reduction = "umap", 
            features = "percent.mt") + 
  scale_colour_gradientn(colours = c("blue","lightblue","yellow","orange","red"))

# UMAP percent.ribo
FeaturePlot(mouse.unannotated,
            reduction = "umap", 
            features = "percent.ribo.protein") + 
  scale_colour_gradientn(colours = c("blue","lightblue","yellow","orange","red"))

# UMAP percent.hb
FeaturePlot(mouse.unannotated,
            reduction = "umap", 
            features = "percent.hb") + 
  scale_colour_gradientn(colours = c("blue","lightblue","yellow","orange","red"))

# UMAP nCount
FeaturePlot(mouse.unannotated,
            reduction = "umap", 
            features = "nCount_RNA") + 
  scale_colour_gradientn(colours = c("blue","lightblue","yellow","orange","red"))

# UMAP nFeature
FeaturePlot(mouse.unannotated,
            reduction = "umap", 
            features = "nFeature_RNA") + 
  scale_colour_gradientn(colours = c("blue","lightblue","yellow","orange","red"))

# UMAP cell.complexity
FeaturePlot(mouse.unannotated,
            reduction = "umap", 
            features = "cell.complexity") + 
  scale_colour_gradientn(colours = c("blue","lightblue","yellow","orange","red"))

# UMAP Ttr expression
FeaturePlot(mouse.unannotated,
            reduction = "umap", 
            features = "Ttr") + 
  scale_colour_gradientn(colours = c("blue","lightblue","yellow","orange","red"))

Violins

VlnPlot(mouse.unannotated,
        features = "nCount_RNA",
        split.by = "seurat_clusters")
## The default behaviour of split.by has changed.
## Separate violin plots are now plotted side-by-side.
## To restore the old behaviour of a single split violin,
## set split.plot = TRUE.
##       
## This message will be shown once per session.

VlnPlot(mouse.unannotated,
        features = "nFeature_RNA",
        split.by = "seurat_clusters")

Cells per cluster

# Cells per sample per cluster
sample_ncells <- FetchData(mouse.unannotated, 
                     vars = c("ident", "sample")) %>%
  dplyr::count(ident,sample) %>%
  tidyr::spread(ident, n)
sample_ncells
##   sample   0    1   2   3   4   5   6   7   8   9  10  11  12  13  14  15  16
## 1   E3.M 580  510 563 263 198 226 291 368 252 288 199 229 123 105 137 108  12
## 2   E3.F 614  372 309 345 306 259 213  31 233 290 126 128 182 146  47  59   8
## 3   E4.M 423  307 552 301 208 289 243 342 245 314 245 187 168 127 105 104 105
## 4   E4.F 820 1027 756 533 419 355 382 320 317 111 347 311 317 205 183 143 283
##    17  18 19 20 21 22  23  24  25
## 1  83 114 43 63 55 54  NA  NA   2
## 2  60  72 86  3 27 44   1  NA   1
## 3  97  66 63 83 59 47   1   1   1
## 4 120  88 90 95 67 51 141 125 120
# Cells per isoform per cluster
isoform_ncells <- FetchData(mouse.unannotated, 
                     vars = c("ident", "isoform")) %>%
  dplyr::count(ident,isoform) %>%
  tidyr::spread(ident, n)
isoform_ncells
##   isoform    0    1    2   3   4   5   6   7   8   9  10  11  12  13  14  15
## 1      E4 1243 1334 1308 834 627 644 625 662 562 425 592 498 485 332 288 247
## 2      E3 1194  882  872 608 504 485 504 399 485 578 325 357 305 251 184 167
##    16  17  18  19  20  21 22  23  24  25
## 1 388 217 154 153 178 126 98 142 126 121
## 2  20 143 186 129  66  82 98   1  NA   3
# Cells per sex per cluster
sex_ncells <- FetchData(mouse.unannotated, 
                     vars = c("ident", "sex")) %>%
  dplyr::count(ident,sex) %>%
  tidyr::spread(ident, n)
sex_ncells
##      sex    0    1    2   3   4   5   6   7   8   9  10  11  12  13  14  15  16
## 1   Male 1003  817 1115 564 406 515 534 710 497 602 444 416 291 232 242 212 117
## 2 Female 1434 1399 1065 878 725 614 595 351 550 401 473 439 499 351 230 202 291
##    17  18  19  20  21  22  23  24  25
## 1 180 180 106 146 114 101   1   1   3
## 2 180 160 176  98  94  95 142 125 121
# Cells per age per cluster
age_ncells <- FetchData(mouse.unannotated, 
                     vars = c("ident", "age")) %>%
  dplyr::count(ident,age) %>%
  tidyr::spread(ident, n)
age_ncells
##         age    0    1    2    3    4    5    6    7    8    9  10  11  12  13
## 1 14 months 2437 2216 2180 1442 1131 1129 1129 1061 1047 1003 917 855 790 583
##    14  15  16  17  18  19  20  21  22  23  24  25
## 1 472 414 408 360 340 282 244 208 196 143 126 124

Gene histogram

# User params
goi <- "Malat1"
threshold <- 1

# Subset data
log2.threshold <- log2(threshold + 0.01)
counts.df <- FetchData(mouse.unannotated, vars = goi)
colnames(counts.df) <- "counts"
log2.counts.df <- log2(counts.df + 0.01)

# Histogram
title <- paste0(goi, "\nnCount_RNA > ", threshold)
hist1 <- ggplot(counts.df, aes(x = counts)) + 
  geom_histogram(bins = 100, fill = "gray", color = "black") + 
  labs(title = title, x=NULL, y=NULL) +
  xlab(paste0(goi, " nCount_RNA")) + ylab("# of Samples") + theme_bw() +
  geom_vline(xintercept = threshold, col = "blue") +
  annotate("rect", 
           xmin = -Inf,
           xmax = threshold, 
           ymin = 0, 
           ymax=Inf, 
           alpha=0.2, 
           fill="chocolate4") +
  annotate("rect", 
           xmin = threshold,
           xmax = Inf, 
           ymin = 0, 
           ymax=Inf, 
           alpha=0.2, 
           fill="deepskyblue")

# Histogram log transformed
hist2 <- ggplot(log2.counts.df, aes(x = counts)) + 
  geom_histogram(bins = 100, fill = "gray", color = "black") + 
  labs(title = title, x=NULL, y=NULL) +
  xlab(paste0("Log2(",goi, " nCount_RNA)")) + ylab("# of Samples") + theme_bw() +
  geom_vline(xintercept = log2.threshold, col = "blue") +
  annotate("rect", 
           xmin = -Inf,
           xmax = log2.threshold, 
           ymin = 0, 
           ymax=Inf, 
           alpha=0.2, 
           fill="chocolate4") +
  annotate("rect", 
           xmin = log2.threshold,
           xmax = Inf, 
           ymin = 0, 
           ymax=Inf, 
           alpha=0.2, 
           fill="deepskyblue")

# plot
plots1 <- list(hist1,hist2)
layout1 <- rbind(c(1),c(2))
grid1 <- grid.arrange(grobs = plots1, layout_matrix = layout1)

Percent gene

# user define variable
goi <- "Malat1"

# Extract counts data
DefaultAssay(mouse.unannotated) <- "RNA"
Idents(mouse.unannotated) <- "SCT_snn_res.0.5"
geneData <- FetchData(mouse.unannotated,
                      vars = goi,
                      slot = "counts")
geneData <- geneData > 1
table(geneData)
## geneData
## FALSE  TRUE 
##   693 20544
mouse.unannotated$Expression <- geneData

FetchData(mouse.unannotated,
          vars = c("ident", "Expression")) %>%
  dplyr::count(ident, Expression) %>%
  tidyr::spread(ident, n)
##   Expression    0    1    2    3    4   5    6    7    8    9  10  11  12  13
## 1      FALSE   NA  147    1    1   NA 273   26   NA   NA   NA  NA  96  NA 112
## 2       TRUE 2437 2069 2179 1441 1131 856 1103 1061 1047 1003 917 759 790 471
##    14  15  16  17  18  19  20  21  22  23  24  25
## 1  19  NA   9   4  NA   5  NA  NA  NA  NA  NA  NA
## 2 453 414 399 356 340 277 244 208 196 143 126 124
# Plot
mouse.unannotated@meta.data %>%
  group_by(SCT_snn_res.0.5, Expression) %>%
  dplyr::count() %>%
  group_by(SCT_snn_res.0.5) %>%
  dplyr::mutate(percent = 100*n/sum(n)) %>%
  ungroup() %>%
  ggplot(aes(x=SCT_snn_res.0.5,y=percent, fill=Expression)) +
  geom_col() +
  ggtitle(paste0("Percentage of cells with > 1 counts for ", goi)) +
  theme(axis.text.x = element_text(angle = 90))

Cluster tree

  • Cluster trees are helpful in deciding what clusters to merge.
mouse.unannotated <- BuildClusterTree(object = mouse.unannotated,
                                     dims = 1:15, # min.pc in processing script
                                     reorder = FALSE,
                                     reorder.numeric = FALSE)
tree <- mouse.unannotated@tools$BuildClusterTree
tree$tip.label <- paste0(tree$tip.label)

ggtree::ggtree(tree, aes(x, y)) +
  scale_y_reverse() +
  ggtree::geom_tree() +
  ggtree::theme_tree() +
  ggtree::geom_tiplab(offset = 1) +
  ggtree::geom_tippoint(color = cluster_colors[1:length(tree$tip.label)], shape = 16, size = 5) +
  coord_cartesian(clip = 'off') +
  theme(plot.margin = unit(c(0,2.5,0,0), 'cm'))

Potential Markers

MACS markers

  • Ptprc - located on most haematopoietic cells, positive selection of leukocytes
  • Pecam1 - adhesion and signaling receptor that is expressed on endothelial and hematopoietic cells, positive selection of endothelial cells
  • Lyve1 - found primarily on lymphatic endothelial cells
VlnPlot(mouse.unannotated,
        features = "Ptprc",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Lyve1",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Pecam1",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

B-cells / Plasma cells

  • Cd19: expressed in B-cells and follicular dendritic cells
  • Fcrla: a B-cell specific protein in mice
  • Cd79a & Cd79b: together form BCR complex
  • Sdc1: Plasma cells
VlnPlot(mouse.unannotated,
        features = "Cd19",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Fcrla",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Cd79a",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Cd79b",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Sdc1",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

T-cells

  • Trac/Cd3d/Cd3e/Cd3g are components of the TCR
VlnPlot(mouse.unannotated,
        features = "Trac",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Cd3e",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Cd8a",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Cd4",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

Endothelial cells

VlnPlot(mouse.unannotated,
        features = "Ly6c1",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Flt1",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

Fibroblasts

VlnPlot(mouse.unannotated,
        features = "Col1a1",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Col1a2",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

Mast cells

VlnPlot(mouse.unannotated,
        features = "Fcer1a",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Kit", # aka Cd117
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

Macrophage

VlnPlot(mouse.unannotated,
        features = "C1qa",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "C1qb",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Mki67",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

Monocytes/DCs

  • Cd11c/Itgax: Monocytes & DCs
VlnPlot(mouse.unannotated,
        features = "Itgax", # Cd11c
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Cd209a",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

Neutrophils

VlnPlot(mouse.unannotated,
        features = "Ly6g",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Retnlg",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

Pericytes & SMCs

VlnPlot(mouse.unannotated,
        features = "Acta2",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Myl9",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Rgs5",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

Schwann cells

VlnPlot(mouse.unannotated,
        features = "Cdh19",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

VlnPlot(mouse.unannotated,
        features = "Mpz",
        cols = cluster_colors,
        split.by = "SCT_snn_res.0.5")

Sandro’s markers

Idents(mouse.unannotated) <- "SCT_snn_res.0.5"
goi1 <- c("Cd3e","Trbc1","Cd4","Cd8a","Foxp3","Tbx21","Gata3","Thy1",
         "Cd19","Ms4a1","Cd27","Ighg1","Ptprc","Ly6g","Itgam","Aif1","Klrb1")
goi2 <- c("Adgre1","Ms4a3","Ly6c2","Mrc1","Lyz2","Cd74","Cd83","Cd14","H2-Aa",
          "H2-Ab1","Sirpa","Xcr1","Siglech","Itgax","Kit","Mcpt4")
goi3 <- c("Pdgfra","Col1a1","Lum","Pdgfrb","Rgs5","Cspg4","Acta2","Tagln",
          "Pecam1","Cd34","Plvap","Stmn2","Slc38a5","Vwf","Mfsd2a","Cldn5")
goi4 <- c("Prox1","Flt4","Pdpn","Lyve1","Sox10","Mbp","Fgf13","Kcnab2","Tubb3",
          "Slc17a6","Shank2","Erbb4","Park7","Kif5b","Slc4a1","Hmbs")
goi5 <- c("Pecam1","Flt4","Itgam","Mrc1","Cd3e","Gata3","Cd19","Ly6g","Pdgfrb",
          "Kit","Col1a1","Pmp22","Hmbs")

v1 <- VlnPlot(mouse.unannotated,
              features = goi1,
              cols = cluster_colors,
              split.by = "SCT_snn_res.0.5",
              flip = TRUE,
              stack = TRUE)
## Warning in FetchData.Seurat(object = object, vars = features, slot = slot): The
## following requested variables were not found: Ighg1
v1

v2 <- VlnPlot(mouse.unannotated,
              features = goi2,
              cols = cluster_colors,
              split.by = "SCT_snn_res.0.5",
              flip = TRUE,
              stack = TRUE)
v2

v3 <- VlnPlot(mouse.unannotated,
              features = goi3,
              cols = cluster_colors,
              split.by = "SCT_snn_res.0.5",
              flip = TRUE,
              stack = TRUE)
v3

v4 <- VlnPlot(mouse.unannotated,
              features = goi4,
              cols = cluster_colors,
              split.by = "SCT_snn_res.0.5",
              flip = TRUE,
              stack = TRUE)
v4

v5 <- VlnPlot(mouse.unannotated,
              features = goi5,
              cols = cluster_colors,
              split.by = "SCT_snn_res.0.5",
              flip = TRUE,
              stack = TRUE)
v5

Automatically detect markers

# Find markers for each cluster
# Compares single cluster vs all other clusters
# Default is logfc.threshold = 0.25, min.pct = 0.5
Idents(mouse.unannotated) <- "seurat_clusters"
all.markers <- FindAllMarkers(object = mouse.unannotated,
                              assay = "RNA",
                              test.use = "MAST",
                              verbose = TRUE)

# add column
all.markers$delta_pct <- abs(all.markers$pct.1 - all.markers$pct.2)

# rename columns and rows
rownames(all.markers) <- 1:nrow(all.markers)
all.markers <- all.markers[,c(6,7,1,5,2:4,8)]
colnames(all.markers)[c(6,7)] <- c("pct_1","pct_2")

# save
saveRDS(all.markers, "../../rObjects/pass1_unannotated_cluster_markers.rds")
# read object
all.markers <- readRDS("../../rObjects/pass1_unannotated_cluster_markers.rds")

# more stringent filtering
all.markers <- all.markers[all.markers$p_val_adj < 0.01,]

# compare 
table(all.markers$cluster)
## 
##    0    1    2    3    4    5    6    7    8    9   10   11   12   13   14   15 
## 2803 2724    0 3565 2517 4124 2598 3533 2796 4615 3234 3100 3558 4186 2745 2520 
##   16   17   18   19   20   21   22   23   24   25 
##  905 2522 2261 2215 3346 2153 1968 1609  513  770
# subset
cluster0 <- all.markers[all.markers$cluster == 0,]
cluster1 <- all.markers[all.markers$cluster == 1,]
cluster2 <- all.markers[all.markers$cluster == 2,]
cluster3 <- all.markers[all.markers$cluster == 3,]
cluster4 <- all.markers[all.markers$cluster == 4,]
cluster5 <- all.markers[all.markers$cluster == 5,]
cluster6 <- all.markers[all.markers$cluster == 6,]
cluster7 <- all.markers[all.markers$cluster == 7,]
cluster8 <- all.markers[all.markers$cluster == 8,]
cluster9 <- all.markers[all.markers$cluster == 9,]
cluster10 <- all.markers[all.markers$cluster == 10,]
cluster11 <- all.markers[all.markers$cluster == 11,]
cluster12 <- all.markers[all.markers$cluster == 12,]
cluster13 <- all.markers[all.markers$cluster == 13,]
cluster14 <- all.markers[all.markers$cluster == 14,]
cluster15 <- all.markers[all.markers$cluster == 15,]
cluster16 <- all.markers[all.markers$cluster == 16,]
cluster17 <- all.markers[all.markers$cluster == 17,]
cluster18 <- all.markers[all.markers$cluster == 18,]
cluster19 <- all.markers[all.markers$cluster == 19,]
cluster20 <- all.markers[all.markers$cluster == 20,]
cluster21 <- all.markers[all.markers$cluster == 21,]
cluster22 <- all.markers[all.markers$cluster == 22,]
cluster23 <- all.markers[all.markers$cluster == 23,]
cluster24 <- all.markers[all.markers$cluster == 24,]
cluster25 <- all.markers[all.markers$cluster == 25,]

# save
write.table(all.markers,
            "../../results/pass_1/all_clusters/markers/unannotated_cluster_markers_pvaladj_0.01.tsv", 
            sep = "\t", quote = FALSE, row.names = FALSE)

Cluster Annotations

Cluster 0: Macrophages 1

VlnPlot(mouse.unannotated,
        features = cluster0$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 1: Ambient RNAs 1

VlnPlot(mouse.unannotated,
        features = cluster1$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 2: Endothelial 1

# no markers
# cluster tree groups with cluster 6
VlnPlot(mouse.unannotated,
        features = "Flt1",
        split.by = "seurat_clusters",
        cols = cluster_colors)

Cluster 3: Fibroblasts 1

VlnPlot(mouse.unannotated,
        features = cluster3$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 4: Macrophages 2

VlnPlot(mouse.unannotated,
        features = cluster4$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 5: Macrophages 3

Cluster 6: Endothelial 2

VlnPlot(mouse.unannotated,
        features = cluster6$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 7: B cells 1

VlnPlot(mouse.unannotated,
        features = cluster7$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 8: Monocytes/DCs

VlnPlot(mouse.unannotated,
        features = cluster8$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 9: Neutrophils

VlnPlot(mouse.unannotated,
        features = cluster9$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 10: T cells 1

VlnPlot(mouse.unannotated,
        features = cluster10$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 11: Endothelial 3

VlnPlot(mouse.unannotated,
        features = cluster11$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 12: Innate lymphoid cells

VlnPlot(mouse.unannotated,
        features = cluster12$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 13: Fibroblasts 2

VlnPlot(mouse.unannotated,
        features = cluster13$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 14: Pericytes & Smooth muscle cells

VlnPlot(mouse.unannotated,
        features = cluster14$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 15: T cells 2

VlnPlot(mouse.unannotated,
        features = cluster15$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 16: Macrophages 4

VlnPlot(mouse.unannotated,
        features = cluster16$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 17: Vwf high endothelial 2

VlnPlot(mouse.unannotated,
        features = cluster17$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 18: Mast cells

VlnPlot(mouse.unannotated,
        features = cluster18$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 19: Macrophages 5

VlnPlot(mouse.unannotated,
        features = cluster19$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 20: B cells 2

VlnPlot(mouse.unannotated,
        features = cluster20$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 21: Schwann cells 1

VlnPlot(mouse.unannotated,
        features = cluster21$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 22: Plasma cells

VlnPlot(mouse.unannotated,
        features = cluster22$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 23: Ambient RNAs 2

  • Pericytes & SMCs (Acta2,Rgs5) and Schwann cells (Mpz,Cdh19)
VlnPlot(mouse.unannotated,
        features = cluster23$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 24: Ambient RNAs 3

VlnPlot(mouse.unannotated,
        features = cluster24$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Cluster 25: Ambient RNAs 4

VlnPlot(mouse.unannotated,
        features = cluster25$gene[1:20],
        cols = cluster_colors,
        stack = TRUE,
        flip = TRUE,
        split.by = "seurat_clusters")

Merge cluster names

# Rename all identities
mouse.annotated <- RenameIdents(object = mouse.unannotated,
                                "0" = "Macrophages",       # Mrc1
                                "1" = "Ambient RNAs",      # Ttr,Grik1,Igha
                                "2" = "Endothelial",       # Flt1
                                "3" = "Fibroblasts",       # Col1a1
                                "4" = "Macrophages",       # C1qb
                                "5" = "Macrophages",       # C1qc
                                "6" = "Endothelial",       # Vwf
                                "7" = "B cells",           # Cd19
                                "8" = "Monocytes and Dendritic cells",    # Fgr
                                "9" = "Neutrophils",       # Retnlg
                                "10" = "T cells",          # Cd3e
                                "11" = "Endothelial",      # Ly6c1
                                "12" = "Innate lymphoid cells",           # Il7r
                                "13" = "Fibroblasts",      # Col1a2
                                "14" = "Pericytes & Smooth muscle cells", # Myh11
                                "15" = "T cells",          # Cd3e
                                "16" = "Macrophages",
                                "17" = "Endothelial",      # Vwf
                                "18" = "Mast cells",       # Kit
                                "19" = "Macrophages",      # Mki67
                                "20" = "B cells",          # Cd19
                                "21" = "Schwann cells",    # Cdh19,Mpz
                                "22" = "Plasma cells",     # Sdc1
                                "23" = "Ambient RNAs",     # Mpz,Sch19,Rgs5,Acta2
                                "24" = "Ambient RNAs",     # Ttr,Enpp2,A2m,Rgs5
                                "25" = "Ambient RNAs")     # Rgs5,Mgp,Flt1
mouse.annotated$merged_clusters <- Idents(mouse.annotated)
# set colors
merged_colors <- c("darkred","firebrick1","gold","darkolivegreen2","darkgreen",
                   "cyan","cornflowerblue","blue","darkorchid1","deeppink",
                   "plum1","burlywood3","azure4", "black")

# umap
umap <- DimPlot(object = mouse.annotated, 
        reduction = "umap",
        repel = TRUE,
        group.by = "merged_clusters",
        cols = merged_colors)
umap

Markers

v <- VlnPlot(mouse.annotated,
        features = c("C1qb","Flt1","Col1a2","Cd19","Itgax","Fgr","Retnlg","Cd3e",
                     "Il7r","Rgs5","Kit","Mpz","Sdc1"),
        split.by = "merged_clusters",
        flip = TRUE,
        stack = TRUE,
        cols = merged_colors)
v

v
path <- "../../results/pass_1/all_clusters/markers/cluster_markers_violin"
saveToPDF(paste0(path, ".pdf"), width = 8, height = 6)
dev.off()

Shiny App

Cleanup object

mouse.annotated@assays$RNA@var.features <- 
  mouse.annotated@assays$SCT@var.features
metadata <- mouse.annotated@meta.data
metadata <- metadata[,c(25,26,2:14)]
colnames(metadata)[1] <- "unannotated_clusters"
mouse.annotated@meta.data <- metadata
mouse.annotated@assays$SCT@meta.features <- metadata
mouse.annotated@assays$RNA@meta.features <- metadata

Output directory

# make shiny folder
setwd("../../results/pass_1/all_clusters")
DefaultAssay(mouse.annotated) <- "RNA"
Idents(mouse.annotated) <- mouse.annotated$merged_clusters
sc.config <- createConfig(mouse.annotated)
makeShinyApp(mouse.annotated, sc.config, gene.mapping = TRUE,
             shiny.title = "Pass 1: Mouse Meningeal Dura scRNAseq")