3C. Annotate the initially identified cell clusters

Export cluster markers

Idents(obj) <- obj$harmony_clusters
features  <- Features(obj)
no_ribo   <- features[!grepl("^Rp[sl][[:digit:]]|^Rplp[[:digit:]]|^Rpsa", features)]
no_ribo_mito <- sort(no_ribo[!grepl("^mt-(.*)", no_ribo)])

all_pos_markers <- FindAllMarkers(obj, only.pos = TRUE, random.seed=99,
                                  test.use="wilcox", features=no_ribo_mito,
                                  min.pct = 0.01, logfc.threshold = 1.0)
all_pos_markers <- Add_Pct_Diff(marker_dataframe = all_pos_markers)

top25_pos <- data.table(all_pos_markers %>% group_by(cluster) %>%
    slice_min(n = 25, order_by = p_val_adj))

fwrite(top25_pos,
    paste0(sample_name, "_Top25PositiveMarkers_", Sys.Date(), ".csv"),
    sep=",")

Annotate cell populations

cluster cell_type assigned_celltype
0 SMC transformed juxtaglomerular cell tSMC_1
1 pericytes PC
2 SMC transformed juxtaglomerular cell tSMC_2
3 SMC transformed juxtaglomerular cell tSMC_3
4 juxtaglomerular cells JG_1
5 proximal tubule epithelial cells PT
6 B cell B
7 proximal tubule segment 3 epithelial cell PTs3
8 juxtaglomerular cells JG_2
9 proximal straight tubule cells PST
10 glomerular endothelial cell gEnC
11 proximal tubule segment 3 epithelial cell PTs3
12 proximal straight tubule cells PST
13 PC transformed tPC
14 Ren1+ parietal epithelial cells PEC
15 unknown unknown
16 macrophage Mp
17 epithelial cell EC
18 transformed mesangial cells tMC
19 transformed mesangial cells tMC
20 intercalated cell IC
21 collecting duct principal cell CDPC
22 proximal tubule segment 3 epithelial cell PTs3
23 T cell T
24 fibroblasts FB
25 proximal tubule segment 3 epithelial cell PTs3
26 proximal tubule epithelial cells PT
27 microphage like smooth muscle cells MpSMC
cluster_identities <- c("tSMC_1", "PC", "tSMC_2", "tSMC_3", "JG_1", "PT", "B",
                        "PTs3", "JG_2", "PST", "gEnC", "PTs3", "PST", "tPC",
                        "PEC", "unknown", "Mp", "EC", "tMC", "tMC", "IC",
                        "CDPC", "PTs3", "T", "FB", "PTs3", "PT", "MpSMC")
names(cluster_identities) <- c("0", "1", "2", "3", "4", "5", "6", "7", "8",
                               "9", "10", "11", "12", "13", "14", "15", "16",
                               "17", "18", "19", "20", "21", "22", "23", "24",
                               "25", "26", "27")

ren1c_named <- copy(obj)
ren1c_named <- RenameIdents(ren1c_named, cluster_identities)
ren1c_named[["assigned_celltype"]] <- Idents(ren1c_named)

pal_named         <- get_random_grid_colors(length(unique(ren1c_named$assigned_celltype)))
names(pal_named)  <- c(unique(ren1c_named$assigned_celltype))

pal_named[["WT"]]  <- uva_blue
pal_named[["KO"]]  <- uva_orange

p1 <- DimPlot(ren1c_named, reduction = "umap_harmony", cols=pal_named,
              group.by = "condition") & 
        custom_theme() &
        theme(legend.position = "right")
p2 <- DimPlot(ren1c_named, reduction = "umap_harmony", cols=pal_named,
              group.by = "assigned_celltype", label = TRUE, repel = TRUE) & 
        custom_theme() &
        theme(legend.position = "right")

ggsave(paste0("figures/Annotated_DimPlot_", Sys.Date(), ".svg"),
        plot = p1 | p2, width = (540*0.0369), height = (180*0.0369),
        units = "mm", create.dir=TRUE)

Visualize renin lineage populations

jg_markers  <- c("Ren1", "Akr1b7")
smc_markers <- c("Acta2", "Smtn", "Tagln", "Myh11", "Cnn1")
pc_markers  <- c("Pdgfrb", "Acta2", "Des", "Cspg4", "Rgs5")
mc_markers  <- c("Serpine2", "Map3k20", "Pde3a", "Hopx", "Myl9", "Filip1l")

plots <- VlnPlot(ren1c_named,
    features = c(jg_markers, smc_markers),
    assay = "RNA", col=pal_named,
    group.by = "assigned_celltype", pt.size = 0, combine = FALSE)

wrap_plots(plots = plots, ncol = 3) & NoLegend()

renin_lineage_cells <- rownames(ren1c_named@meta.data[ren1c_named@meta.data$assigned_celltype %in% c("tSMC_1", "tSMC_2", "tSMC_3", "PC", "tPC", "JG_1", "JG_2", "MpSMC"),])
renin_lineage <- subset(ren1c_named, cells=renin_lineage_cells)
renin_lineage$assigned_celltype <- factor(renin_lineage$assigned_celltype, 
    levels = c("JG_1", "JG_2", "tSMC_1", "tSMC_2", "tSMC_3", "MpSMC", "PC", "tPC"))

stk_vln_renin <- Stacked_VlnPlot(renin_lineage,
    features = c(jg_markers, smc_markers),
    assay = "RNA", colors_use=pal_named,
    group.by = "assigned_celltype", pt.size = 0, combine = T)

stk_vln_renin_split <- Stacked_VlnPlot(renin_lineage,
    features = c(jg_markers, smc_markers),
    assay = "RNA", 
    split.by = "condition", plot_legend=T,
    group.by = "assigned_celltype", pt.size = 0, combine = T)

ggsave(paste0("figures/ReninLineageCells_JG-SMC-Markers_StackedVln_", Sys.Date(), ".svg"),
        plot = stk_vln_renin | stk_vln_renin_split,
        width = (540*0.0369), height = (540*0.0369), units="mm")

Identify JGs and SMCs

gxm_mat <- LayerData(ren1c_named, assay="RNA") # OR: ren1c_named@assays$RNA$data
if(exists("jg_names")){rm(jg_names)}
for (gene in jg_markers) {
    if (gene == "Ren1") {jg_cutoff <- 0.95}
    if (gene == "Akr1b7") {jg_cutoff <- 0.90}
    jg_gxm <- gxm_mat[grepl(paste0("^", gene, "$"), rownames(gxm_mat)),]
    jg_upper_quartile <- names(jg_gxm[jg_gxm > quantile(gxm_mat[grep(paste0("^", gene, "$"), rownames(gxm_mat)),], jg_cutoff)[[1]]])
    message(paste0(gene, ": ", length(jg_upper_quartile)))
    if (exists("jg_names")) {
        jg_names <- Reduce(intersect, list(jg_names, jg_upper_quartile))
    } else {
        jg_names <- jg_upper_quartile
    }
}
JG_index  <- grep(paste(jg_names, collapse = "|"), names(Idents(ren1c_named)))

if(exists("smc_names")){rm(smc_names)}
for (gene in c(smc_markers, jg_markers)) {
    if (gene == "Acta2") {smc_cutoff <- 0.75}
    if (gene == "Smtn") {smc_cutoff <- 0.75}
    if (gene == "Tagln") {smc_cutoff <- 0.50}
    if (gene == "Myh11") {smc_cutoff <- 0.50}
    if (gene == "Cnn1") {smc_cutoff <- 0.50}
    if (gene == "Ren1") {smc_cutoff <- 0.05}
    if (gene == "Akr1b7") {smc_cutoff <- 0.05}
    smc_gxm <- gxm_mat[grepl(paste0("^", gene, "$"), rownames(gxm_mat)),]
    if (gene == "Ren1" | gene == "Akr1b7") {
        smc_upper_quartile <- names(smc_gxm[smc_gxm <= quantile(gxm_mat[grep(paste0("^", gene, "$"), rownames(gxm_mat)),], smc_cutoff)[[1]]])
    } else {
        smc_upper_quartile <- names(smc_gxm[smc_gxm > quantile(gxm_mat[grep(paste0("^", gene, "$"), rownames(gxm_mat)),], smc_cutoff)[[1]]])
    }
    message(paste0(gene, ": ", length(smc_upper_quartile)))
    if (exists("smc_names")) {
        if (gene == "Ren1" | gene == "Akr1b7") {
            #smc_names <- smc_names[!(smc_names %in% smc_upper_quartile)]
            smc_names <- Reduce(intersect, list(smc_names, smc_upper_quartile))
        } else {
            smc_names <- Reduce(intersect, list(smc_names, smc_upper_quartile))
        }
    } else {
        smc_names <- smc_upper_quartile
    }
}

SMC_names <- rownames(ren1c_named@meta.data[ren1c_named@meta.data$assigned_celltype %in% c("tSMC_1", "tSMC_2", "tSMC_3"),])
gene <- "Ren1"
ren1_gxm <- gxm_mat[grepl(paste0("^", gene, "$"), rownames(gxm_mat)),]
SMC_ren1_gxm <- ren1_gxm[SMC_names]
no_renin_smc <- names(SMC_ren1_gxm[SMC_ren1_gxm==0])

gene <- "Akr1b7"
akr1b7_gxm <- gxm_mat[grepl(paste0("^", gene, "$"), rownames(gxm_mat)),]
SMC_akr1b7_gxm <- akr1b7_gxm[SMC_names]
quantile(akr1b7_gxm)
no_akr1b7_smc <- names(SMC_akr1b7_gxm[SMC_akr1b7_gxm==0])

true_smc_names <- intersect(no_renin_smc, no_akr1b7_smc)
SMC_cells <- list(SMCs = true_smc_names)

Idents(ren1c_named) <- ren1c_named$assigned_celltype 
SMC_index  <- grep(paste(true_smc_names, collapse = "|"), names(Idents(ren1c_named)))
ren1c_named <- SetIdent(ren1c_named, cells = SMC_index, value = 'SMC')
ren1c_named[["assigned_celltype"]] <- Idents(ren1c_named)

pal_named[["SMC"]] <- uva_green

ggsave(paste0("ren1c-1mo_annotated-UMAP-wSMC_", Sys.Date(), ".svg"),
       plot = DimPlot_scCustom(ren1c_named, reduction = "umap_harmony",
                colors_use=pal_named, group.by = "assigned_celltype",
                label = TRUE, repel = TRUE) & 
              custom_theme() &
              theme(legend.position = "right"))

renin_lineage_cells <- rownames(ren1c_named@meta.data[ren1c_named@meta.data$assigned_celltype %in% c("SMC", "tSMC_1", "tSMC_2", "tSMC_3", "MpSMC", "PC", "tPC", "JG_1", "JG_2"),])

JG_names <- rownames(ren1c_named@meta.data[ren1c_named@meta.data$assigned_celltype %in% c("JG_1", "JG_2"),])

gene <- "Acta2"
acta2_gxm <- gxm_mat[grepl(paste0("^", gene, "$"), rownames(gxm_mat)),]
acta2_median_exp <- median(acta2_gxm[renin_lineage_cells])
JG_acta2_gxm <- acta2_gxm[JG_names]
low_acta2_JG <- names(JG_acta2_gxm[JG_acta2_gxm < acta2_median_exp])

gene <- "Smtn"
smtn_gxm <- gxm_mat[grepl(paste0("^", gene, "$"), rownames(gxm_mat)),]
smtn_median_exp <- median(smtn_gxm[renin_lineage_cells])
JG_smtn_gxm <- smtn_gxm[JG_names]
low_smtn_JG <- names(JG_smtn_gxm[JG_smtn_gxm < smtn_median_exp])

gene <- "Tagln"
tagln_gxm <- gxm_mat[grepl(paste0("^", gene, "$"), rownames(gxm_mat)),]
tagln_median_exp <- median(tagln_gxm[renin_lineage_cells])
JG_tagln_gxm <- tagln_gxm[JG_names]
low_tagln_JG <- names(JG_tagln_gxm[JG_tagln_gxm < tagln_median_exp])

gene <- "Myh11"
myh11_gxm <- gxm_mat[grepl(paste0("^", gene, "$"), rownames(gxm_mat)),]
myh11_median_exp <- median(myh11_gxm[renin_lineage_cells])
JG_myh11_gxm <- myh11_gxm[JG_names]
low_myh11_JG <- names(JG_myh11_gxm[JG_myh11_gxm < myh11_median_exp])

gene <- "Cnn1"
cnn1_gxm <- gxm_mat[grepl(paste0("^", gene, "$"), rownames(gxm_mat)),]
cnn1_median_exp <- median(cnn1_gxm[renin_lineage_cells])
JG_cnn1_gxm <- cnn1_gxm[JG_names]
low_cnn1_JG <- names(JG_cnn1_gxm[JG_cnn1_gxm < cnn1_median_exp])

true_jg_names <- Reduce(intersect, list(low_acta2_JG,
                                        low_smtn_JG,
                                        low_tagln_JG,
                                        low_myh11_JG,
                                        low_cnn1_JG))

JG_cells <- list(JGs = true_jg_names)
Cell_Highlight_Plot(seurat_object = ren1c_named, reduction = "umap_harmony",
    cells_highlight = JG_cells, label=TRUE) +
    custom_theme() +
    theme(legend.position = "right")

Idents(ren1c_named) <- ren1c_named$assigned_celltype 
JG_index  <- grep(paste(true_jg_names, collapse = "|"), names(Idents(ren1c_named)))
ren1c_named <- SetIdent(ren1c_named, cells = JG_index, value = 'JG')
ren1c_named[["assigned_celltype"]] <- Idents(ren1c_named)

pal_named[["JG"]] <- uva_magenta

Idents(ren1c_named) <- ren1c_named$assigned_celltype 
ren1c_named <- RenameIdents(object = ren1c_named, `JG_2` = "JG")
ren1c_named <- RenameIdents(object = ren1c_named, `JG_1` = "tSMC_4")
ren1c_named[["assigned_celltype"]] <- Idents(ren1c_named)
pal_named[["tSMC_4"]] <- uva_teal

ggsave(paste0("ren1c_named_annotated-UMAP-wSMC-wJG_", Sys.Date(), ".svg"),
       plot = DimPlot_scCustom(ren1c_named, reduction = "umap_harmony",
                  colors_use=pal_named, group.by = "assigned_celltype",
                  label = TRUE, repel = TRUE) & 
              custom_theme() &
              theme(legend.position = "right"))
  • tSMC_4 is the most renin-like of the transitional SMCs and becomes 3
  • tSMC_2 and 3 should combine
  • Overall, this structure looks like a transition between mature SMCs to JGs
Idents(ren1c_named) <- ren1c_named$assigned_celltype
ren1c_named <- RenameIdents(object = ren1c_named, `tSMC_3` = "tSMC_2")
ren1c_named <- RenameIdents(object = ren1c_named, `tSMC_4` = "tSMC_3")
ren1c_named[["assigned_celltype"]] <- Idents(ren1c_named)

pal_named[["JG"]] <- uva_magenta
pal_named[["SMC"]] <- uva_green

ggsave(paste0("ren1c_named_final-annotated-UMAP-wSMC-wJG_", Sys.Date(), ".svg"),
        plot = DimPlot_scCustom(ren1c_named, reduction = "umap_harmony",
                    colors_use=pal_named, group.by = "assigned_celltype",
                    label = TRUE, repel = TRUE) & 
                custom_theme() &
                theme(legend.position = "right"))

renin_lineage_cells <- rownames(ren1c_named@meta.data[ren1c_named@meta.data$assigned_celltype %in% c("tSMC_1", "tSMC_2", "tSMC_3", "PC", "tPC", "JG", "SMC", "MpSMC"),])
renin_lineage <- subset(ren1c_named, cells=renin_lineage_cells)
renin_lineage$assigned_celltype <- factor(renin_lineage$assigned_celltype, 
    levels = c("SMC", "tSMC_1", "tSMC_2", "tSMC_3", "MpSMC", "JG", "PC", "tPC"))

stk_vln_renin <- Stacked_VlnPlot(renin_lineage,
    features = c(jg_markers, smc_markers),
    assay = "RNA", colors_use=pal_named,
    group.by = "assigned_celltype", pt.size = 0, combine = T)

stk_vln_renin_split <- Stacked_VlnPlot(renin_lineage,
    features = c(jg_markers, smc_markers),
    assay = "RNA", colors_use=c(uva_orange, uva_blue),
    split.by = "condition", plot_legend=T,
    group.by = "assigned_celltype", pt.size = 0, combine = T)

ggsave(paste0("figures/ReninLineageCells_1mo_JG-SMC-Markers_StackedVln_RenameIdents_", Sys.Date(), ".svg"),
        plot = stk_vln_renin | stk_vln_renin_split,
        width = (540*0.0369), height = (540*0.0369), units="mm")

all_cells_1mo <- data.table(cell_name = names(ren1c_named$condition),
                            sample_source = ren1c_named$condition,
                            cell_type = Idents(ren1c_named))
all_cells_1mo$sample_source <- ordered(all_cells_1mo$sample_source,
                                       levels = c("WT", "KO"))
all_cells_1mo$cell_type <- factor(all_cells_1mo$cell_type,
    levels = c("SMC", "tSMC_1", "tSMC_2", "tSMC_3", "JG", "PC", "tPC",
               "MpSMC", "B", "Mp", "T", "FB", "PT", "PTs3", "PST", "PTbb",
               "PTs2s3", "CDPC", "IC", "gEnC", "EC", "PEC"))

fwrite(as.data.table(table(all_cells_1mo$sample_source)),
       paste0(sample_name, "_cells-by-source.csv"), sep=",")
fwrite(as.data.table(table(all_cells_1mo$cell_type)),
       paste0(sample_name, "_cells-by-annotation.csv"), sep=",")

Subcluster the renin lineage cells

renin_lineage_cells <- rownames(ren1c_named@meta.data[ren1c_named@meta.data$assigned_celltype %in% c("tSMC_1", "tSMC_2", "tSMC_3", "PC", "tPC", "JG", "SMC"),])
renin_lineage <- subset(ren1c_named, cells=renin_lineage_cells)
renin_lineage$assigned_celltype <- factor(renin_lineage$assigned_celltype, 
    levels = c("SMC", "tSMC_1", "tSMC_2", "tSMC_3", "MpSMC", "JG", "PC", "tPC"))

renin_lineage <- FindNeighbors(renin_lineage, reduction = "integrated_cca", dims = 1:25)
renin_lineage <- FindClusters(renin_lineage, resolution = 0.75, cluster.name = "cca_clusters",
                    random.seed = 99, algorithm=1, group.singletons=TRUE)
renin_lineage <- RunUMAP(renin_lineage, reduction = "integrated_cca", dims = 1:25,
               reduction.name = "umap_cca", seed.use = 99, n.neighbors = 15L,
               n.epochs = 500, min.dist=0.01, spread=1, verbose = FALSE,
               local.connectivity = 1.5L, repulsion.strength = 1, negative.sample.rate = 25L)

renin_lineage <- FindNeighbors(renin_lineage, reduction = "integrated_rpca", dims = 1:25)
renin_lineage <- FindClusters(renin_lineage, resolution = 0.5, cluster.name = "rpca_clusters",
                    random.seed = 99, algorithm=2, group.singletons=TRUE)
renin_lineage <- RunUMAP(renin_lineage, reduction = "integrated_rpca", dims = 1:25,
               reduction.name = "umap_rpca", seed.use = 99, n.neighbors = 15L,
               n.epochs = 500, min.dist=0.1, spread=1.0, verbose = FALSE,
               negative.sample.rate = 25L, repulsion.strength = 1)

renin_lineage <- FindNeighbors(renin_lineage, reduction = "harmony", dims = 1:25)
renin_lineage <- FindClusters(renin_lineage, resolution = 0.75, cluster.name = "harmony_clusters",
                    random.seed = 99, algorithm=1, group.singletons=TRUE)
renin_lineage <- RunUMAP(renin_lineage, reduction = "harmony", dims = 1:25,
               reduction.name = "umap_harmony", seed.use = 99,
               n.neighbors = 15L, n.epochs = 500, min.dist=0.1,
               spread=1.0, verbose = FALSE, negative.sample.rate = 15L)

cca_umap <- DimPlot_scCustom(
  renin_lineage,
  reduction = "umap_cca",
  group.by = c("assigned_celltype", "cca_clusters"),
  combine = FALSE, label.size = 2
)
cca_umap[[1]] | cca_umap[[2]]
rpca_umap <- DimPlot(
  renin_lineage,
  reduction = "umap_rpca",
  group.by = c("assigned_celltype", "rpca_clusters"),
  combine = FALSE, label.size = 2
)
rpca_umap[[1]] | rpca_umap[[2]]
harmony_umap <- DimPlot(
  renin_lineage,
  reduction = "umap_harmony",
  group.by = c("assigned_celltype", "harmony_clusters"),
  combine = FALSE, label.size = 2
)
harmony_umap[[1]] | harmony_umap[[2]]

wrap_plots(c(cca_umap, rpca_umap, harmony_umap), ncol = 3, byrow = F)

sample_name <- "renin-lineage-only"

Idents(renin_lineage) <- renin_lineage$cca_clusters
pal_lineage <- get_random_grid_colors(length(unique(renin_lineage$cca_clusters)))
names(pal_lineage)  <- c(unique(renin_lineage$cca_clusters))

pal_lineage[["WT"]]  <- uva_blue
pal_lineage[["KO"]]  <- uva_orange

p1 <- DimPlot_scCustom(renin_lineage, reduction = "umap_cca",
                       group.by = "assigned_celltype") & 
        custom_theme() &
        theme(legend.position = "right")
p2 <- DimPlot_scCustom(renin_lineage, reduction = "umap_cca",
                       colors_use=pal_lineage,
                       group.by = "cca_clusters", label = TRUE, repel = TRUE) & 
        custom_theme() &
        theme(legend.position = "right")

ggsave(paste0("figures/", sample_name, "_CCA-UMAP_", Sys.Date(), ".svg"),
       plot = p1 | p2, units="mm")
saveRDS(renin_lineage, "renin_lineage_object.rds")

Identify markers for annotation

Idents(renin_lineage) <- renin_lineage$cca_clusters
features  <- Features(renin_lineage)
no_ribo   <- features[!grepl("^Rp[sl][[:digit:]]|^Rplp[[:digit:]]|^Rpsa", features)]
no_ribo_mito <- sort(no_ribo[!grepl("^mt-(.*)", no_ribo)])
no_ribo_mito_hemo <- sort(no_ribo_mito[!grepl("^H2-(.*)", no_ribo_mito)])
filtered_features <- sort(no_ribo_mito_hemo[!grepl("^Hb[ab]-(.*)", no_ribo_mito_hemo)])

lineage_markers <- FindAllMarkers(renin_lineage, only.pos = TRUE, random.seed=99,
                                  test.use="MAST", features=filtered_features,
                                  min.pct = 0.10, logfc.threshold = 0.25)
lineage_markers <- Add_Pct_Diff(marker_dataframe = lineage_markers)

top25_pos <- data.table(lineage_markers %>% group_by(cluster) %>%
    slice_min(n = 25, order_by = p_val_adj))
fwrite(top25_pos,
    paste0(sample_name, "_Top25PositiveMarkers_", Sys.Date(), ".csv"),
    sep=",")

renal_fibrosis_markers <- unique(c("Adm", "Adora2a", "Ager", "Agt", "Agtrap",
    "Akt1", "Alb", "Angpt2", "Anxa1", "Apln", "Aplnr", "Apoe", "Azgp1",
    "Bdkrb1", "Bdkrb2", "Bgn", "Bicc1", "Bmp1", "Bmp2", "Bmp3", "Bmp7",
    "Bmpr1a", "C6", "Cast", "Cat", "Cav1", "Ccl2", "Ccl5", "Ccr7", "Cd248",
    "Cdc42", "Cdh16", "Cebpa", "Cftr", "Clu", "Col1a1", "Col3a1", "Col4a3",
    "Creb1", "Csf1", "Ctnnb1", "Cubn", "Cul3", "Cxcl16", "Cybb", "Cys1",
    "Ddah1", "Dnmt1", "Edn1", "Ednra", "Ednrb", "Efnb2", "Egfr", "Egln1", 
    "Eng", "Epas1", "Epo", "Erbb4", "Fan1", "Fermt2", "Fgf23", "Fn1", "Gli1",
    "Gli2", "Glis2", "Havcr1", "Hgf", "Hipk2", "Hmox1", "Hnrnpf", "Ift140",
    "Ift20", "Igf1", "Ilk", "Invs", "Itga1", "Itgb1", "Jak3", "Kcnj1", "Kcnn4",
    "Kl", "Klf5", "Klhl3", "Lcn2", "Lepr", "Lrp2", "Mmp14", "Mmp2", "Mpz",
    "Mrc2", "Mtor", "N4bp2l2", "Ncf1", "Nfkb1", "Nos3", "Notch1", "Notch3",
    "Nox4", "Nr0b2", "Nr1h4", "Pax8", "Pdpn", "Pkd1", "Pkhd1", "Plat", "Plaur",
    "Postn", "Ppara", "Pparg", "Ppargc1a", "Prkaa1", "Prkcb", "Prkg1", "Prss23",
    "Ptgis", "Ptgs2", "Rasal1", "Rbpj", "Rictor", "Rock1", "Rptor", "Rtn1",
    "S100a4", "Sdc4", "Serpine1", "Shh", "Shroom3", "Sirt1", "Smad7", "Smurf1",
    "Smurf2", "Snai1", "Snai2", "Socs1", "Srebf1", "Stk11", "Tet3", "Tgfb1",
    "Tgfb2", "Tgfb3", "Tgm2", "Thbs1", "Thg1l", "Tlr4", "Trpc3", "Tsc1",
    "Twist1", "Ucp2", "Umod", "Wfdc2", "Wnk4", "Wnt4", "Xpnpep1", "Xylt2"))

renin_lineage <- AddModuleScore_UCell(renin_lineage,
    features=list(Fibrosis=renal_fibrosis_markers))
fibrosis_idents_plt <- VlnPlot_scCustom(renin_seurat,
    features = "Fibrosis_UCell", alpha=0.1, group.by = "cca_clusters",
    split.by="genotype", pt.size = 0, color_seed = 99) +
          custom_theme() +
          ggforce::geom_sina(seed=99, binwidth=100, adjust=4) +
          theme(plot.title = element_text(size = 15/.pt),
                axis.text.x = element_text(size = 12/.pt, angle = 90,
                                           hjust = 1, vjust = 0.5),
                axis.text.y = element_text(size = 12/.pt),
                axis.title.y = element_text(size = 15/.pt)) +
          labs(x="", y="Fibrosis Score (normalized U Statistic)")
ggsave(paste0("figures/", "fibrosis-score_", Sys.Date(), ".png"),
       plot = fibrosis_idents_plt, width = 90, height = 90,
       units = "mm", dpi = 300)
ggsave(paste0("figures/", "fibrosis-score_", Sys.Date(), ".svg"),
       plot = fibrosis_idents_plt, width = 90, height = 90,
       units = "mm", dpi = 300)

Assign labels for the renin lineage populations

cluster_identities <- c("SMC", "tSMC_2", "PC", "fbJG", "tSMC_1",
                        "tSMC_3", "JG", "tPC", "MC")
names(cluster_identities) <- c("0", "1", "2", "3", "4", "5", "6", "7", "8")

lineage_named <- copy(renin_lineage)
lineage_named <- RenameIdents(lineage_named, cluster_identities)
lineage_named[["cell_type"]] <- Idents(lineage_named)

pal_lineage[["JG"]]  <- uva_magenta
pal_lineage[["SMC"]] <- uva_green
pal_lineage[["tSMC_1"]] <- uva_blue
pal_lineage[["tSMC_2"]] <- uva_cyan
pal_lineage[["tSMC_3"]] <- uva_yellow
pal_lineage[["fbJG"]] <- uva_light_orange
pal_lineage[["tPC"]] <- uva_teal
pal_lineage[["MC"]]  <- uva_light_blue
pal_lineage[["PC"]]  <- uva_orange

saveRDS(pal_lineage, "pal_lineage.rds")

dimplot_renamed <- DimPlot_scCustom(lineage_named, reduction = "umap_cca",
        colors_use=pal_lineage, group.by = "cell_type",
        label = TRUE, repel = TRUE) & 
    custom_theme() &
    theme(legend.position = "right")

ggsave(paste0("figures/", sample_name, "_label-UMAP_", Sys.Date(), ".png"),
       plot = dimplot_renamed, width = 180, height = 180,
       units = "mm", dpi = 600)
ggsave(paste0("figures/", sample_name, "_label-UMAP_", Sys.Date(), ".svg"),
       plot = dimplot_renamed, width = 180, height = 180,
       units = "mm", dpi = 300) 

saveRDS(lineage_named, "renin_lineage_named_object.rds")

stk_vln_renin <- Stacked_VlnPlot(lineage_named,
    features = c(jg_markers, smc_markers),
    assay = "RNA", colors_use=pal_lineage,
    group.by = "cell_type", pt.size = 0, combine = T)

stk_vln_renin_split <- Stacked_VlnPlot(lineage_named,
    features = c(jg_markers, smc_markers),
    assay = "RNA", colors_use=c(uva_orange, uva_blue),
    split.by = "condition", plot_legend=T,
    group.by = "cell_type", pt.size = 0, combine = T)

ggsave(paste0("figures/", sample_name, "_ReninMarkers-Vln_", Sys.Date(), ".png"),
       plot = stk_vln_renin | stk_vln_renin_split, width = 360, height = 180, 
       units = "mm", dpi = 600)
ggsave(paste0("figures/", sample_name, "_ReninMarkers-Vln_", Sys.Date(), ".svg"),
       plot = stk_vln_renin | stk_vln_renin_split, width = 360, height = 180,
       units = "mm", dpi = 300) 

Plot renin and smooth muscle lineage markers on reclustered populations

  • Related to Supplemental Figure 3.
renin_seurat <- subset(lineage_named,
    idents = c("SMC", "tSMC_1", "tSMC_2", "tSMC_3", "fbJG", "JG"), invert = F)
renin_seurat$genotype <- factor(renin_seurat$genotype, levels=c("KO", "WT"))

renin_lineage_marker_plt <- VlnPlot_scCustom(renin_seurat,
    features = c("Ren1", "Akr1b7", "Acta2", "Smtn", "Tagln", "Myh11"),
    alpha=0.1, group.by = "cell_type", split.by="genotype", 
    pt.size = 0, color_seed = 99) &
    custom_theme() &
    ggforce::geom_sina(seed=99, binwidth=100, adjust=4) &
    theme(plot.title = element_text(size = 15/.pt),
        axis.text.x = element_text(size = 12/.pt, angle = 90,
                                   hjust = 1, vjust = 0.5),
        axis.text.y = element_text(size = 12/.pt),
        axis.title.y = element_text(size = 15/.pt))

ggsave(paste0("figures/", "ren1c_trajectory-only_renin-lineage-markers_", Sys.Date(), ".png"),
       plot = renin_lineage_marker_plt, width = 270, height = 180, units = "mm", dpi = 300)
ggsave(paste0("figures/", "ren1c_trajectory-only_renin-lineage-markers_", Sys.Date(), ".svg"),
       plot = renin_lineage_marker_plt, width = 270, height = 180, units = "mm", dpi = 300)

### Differential expression — JG versus transformed populations

**Produces manuscript Figure 2b** (volcano plot of DEGs between JG cells and the transformed populations tSMC_1, tSMC_2, tSMC_3, and fbJG in the Ren1c model).

`plotVolcano` is defined in `extra_functions.R`. It calls `FindMarkers` internally (MAST test), writes the DEG CSV to `outdir`, and returns an `EnhancedVolcano` plot object with the top up/down genes and any `highlight_genes` labeled.

source("extra_functions.R")

Idents(renin_seurat) <- renin_seurat$cell_type

features <- rownames(renin_seurat) features <- features[!grepl("^Rp[sl][[:digit:]]|^Rplp[[:digit:]]|^Rpsa", features)] features <- features[!grepl("^mt-", features)] features <- features[!grepl("^Hb[a-b]", features)] features <- features[!grepl("^[0-9]", features)]

p_volcano <- plotVolcano( renin_seurat, group.1 = "JG", group.2 = c("tSMC_1", "tSMC_2", "tSMC_3", "fbJG"), test = "MAST", gene_list = features, no_shape = TRUE, outdir = "figures/DEGs", highlight_genes = c("Ren1", "Akr1b7", "Acta2", "Col4a2", "Cdh13") )

if (!is.null(p_volcano)) { ggsave(paste0("figures/", sample_name, "JG-vs-Transformed_volcano", Sys.Date(), ".svg"), plot = p_volcano, width = 180, height = 180, units = "mm", dpi = 300) }