Use Seurat to find specific markers of each cell type (modified from Dan Skelly)

find_specific_markers.R

library(tidyverse)
library(Seurat)
# Find specific markers for each cell type
# Find specific markers for each cell type

celltypes <- unique(expression_seurat_subset@ident)  # list of cell types (unique ident labels of Seurat object)
obj <-   expression_seurat_subset      # Seurat object
specific_markers <- NULL

# First we do all pairwise comparisons and retain any markers that 
# are even somewhat higher in the cell type of interest
for (celltype1 in celltypes) {
  print("Now working on cluster...")
  print(celltype1)
  other_types <- setdiff(celltypes, celltype1)
  celltype_markers <- NULL
  for (celltype2 in other_types) {
    markers <- FindMarkers(obj, ident.1=celltype1, ident.2=celltype2,
                           only.pos=TRUE,
                           test.use="roc") %>%
      rownames_to_column("gene") %>%
      mutate(ident.2=celltype2)
    celltype_markers <- rbind(celltype_markers, markers)
    markers
  }
  celltype_markers$ident.1 <- celltype1
  specific_markers <- rbind(specific_markers, celltype_markers)
}

write.table(specific_markers,
            file="../Marker/Specific_markers_custom.tsv",
            sep="\t",
            row.names=FALSE,
            col.names=TRUE,
            quote=FALSE)

# In specific_markers data.frame, ident.1 is which cell pop the marker is tagging
# and ident.2 is which other cell pop we are contrasting with
#
# We want markers that are considerably higher in ident.1 than in any other pop.
# As a basic filter we require any marker for ident.1 to be expressed in <50% of 
# cells of each of the other cell types.
# 
# Here we will filter to get the top 10 markers for each cell type

topN <- 5
npop <- length(celltypes)
final_markers <- mutate(specific_markers, pct_diff=pct.1 - pct.2) %>%
    filter(n() == npop - 1, all(pct.2 < 0.5)) %>%
    summarize(mean_AUC=mean(myAUC)) %>%
    mutate(cluster=ident.1) %>%
    filter(mean_AUC > 0.65) %>%
    group_by(cluster) %>%
    arrange(desc(mean_AUC)) %>%
    do(head(., topN))

write.table(specific_markers,
            file="../Marker/Final_markers_custom.tsv",
            sep="\t",
            row.names=FALSE,
            col.names=TRUE,
            quote=FALSE)

png(filename="../Marker/Dotplot_custom_markers_top5.png", width=1600, height=1400, bg = "white", res = 150)
DotPlot(expression_seurat_subset,
        unique(final_markers$gene),
        x.lab.rot = TRUE)
dev.off()