网站仿站教程,网站改了标题会怎么样,网站左侧悬浮导航代码,深圳福田网站建设公司文章目录 介绍数据和代码图1图2图6附图2附图3附图4附图5附图6 介绍
文章提供画图代码和数据#xff0c;本文记录
数据和代码
数据可从以下链接下载#xff08;画图所需要的所有数据#xff09;#xff1a; 百度云盘链接: https://pan.baidu.com/s/1peU1f8_TG2kUKXftkpYq… 文章目录 介绍数据和代码图1图2图6附图2附图3附图4附图5附图6 介绍
文章提供画图代码和数据本文记录
数据和代码
数据可从以下链接下载画图所需要的所有数据 百度云盘链接: https://pan.baidu.com/s/1peU1f8_TG2kUKXftkpYqug 提取码: 7pjy
图1 #### Figure 1: Census of cell types of the mouse uterine tube ######## Packages Load ####library(dplyr)
library(patchwork)
library(Seurat)
library(harmony)
library(ggplot2)
library(cowplot)
library(SoupX)
library(DoubletFinder)
library(data.table)
library(parallel)
library(tidyverse)
library(SoupX)
library(ggrepel)library(ggplot2)
library(gplots)
library(RColorBrewer)
library(viridisLite)
library(Polychrome)
library(circlize)
library(NatParksPalettes)#### Distal and Proximal Datasets ####Distal - readRDS(file ../dataset/Distal_Filtered_Cells.rds , refhook NULL)Proximal - readRDS( file ../dataset/Proximal_Filtered_Cells.rds , refhook NULL)#### Figure 1b: Cells of the Distal Uterine Tube ####Distal_Named - RenameIdents(Distal, 0 Fibroblast 1, 1 Fibroblast 2, 2 Secretory Epithelial,3 Smooth Muscle, 4 Ciliated Epithelial 1, 5 Fibroblast 3, 6 Stem-like Epithelial 1,7 Stem-like Epithelial 2,8 Ciliated Epithelial 2, 9 Blood Endothelial, 10 Pericyte, 11 Intermediate Epithelial, 12 T/NK Cell, 13 Epithelial/Fibroblast, 14 Macrophage, 15 Erythrocyte, 16 Luteal,17 Mesothelial,18 Lymph Endothelial/Epithelial) # Remove cluster due few data points and suspected doubletDistal_Namedactive.ident - factor(x Distal_Namedactive.ident, levels c(Fibroblast 1,Fibroblast 2,Fibroblast 3,Smooth Muscle,Pericyte,Blood Endothelial,Lymph Endothelial/Epithelial,Epithelial/Fibroblast,Stem-like Epithelial 1,Stem-like Epithelial 2,Intermediate Epithelial,Secretory Epithelial,Ciliated Epithelial 1,Ciliated Epithelial 2,T/NK Cell,Macrophage,Erythrocyte,Mesothelial,Luteal))Distal_Named - SetIdent(Distal_Named, value Distal_Namedactive.ident)Fibroblasts - c(#FF9D00 , #FFB653 , #FFCB9A) # Oranges
Muscle - c(#E55451 , #FFB7B2) # Reds
Endothelial - c(#A0E6FF) # Reds
FiboEpi - #FFE0B3 # Reddish Brown
Epi -c(#6E3E6E,#8A2BE2,#CCCCFF,#DA70D6,#DF73FF,#604791) # Blues/Purples
Immune - c( #5A5E6B , #B8C2CC , #FC86AA) # Yellowish Brown
Meso - #9EFFFF # Pink
Lut - #9DCC00 # Greencolors - c(Fibroblasts, Muscle, Endothelial, FiboEpi, Epi, Immune, Meso, Lut)pie(rep(1,length(colors)), colcolors) Distal_Named - subset(Distal_Named, idents c(Fibroblast 1,Fibroblast 2,Fibroblast 3,Smooth Muscle,Pericyte,Blood Endothelial,Epithelial/Fibroblast,Stem-like Epithelial 1,Stem-like Epithelial 2,Intermediate Epithelial,Secretory Epithelial,Ciliated Epithelial 1,Ciliated Epithelial 2,T/NK Cell,Macrophage,Erythrocyte,Mesothelial,Luteal))p1 - DimPlot(Distal_Named,reductionumap,colscolors,pt.size 0.5,label.size 4,label.color black,repel TRUE,labelF) NoLegend() labs(xUMAP_1,yUMAP_2)LabelClusters(p1, idident, color black, repel T , size 4, box.padding .75)ggsave(filename FIG1b_all_distal_umap.pdf, plot p1, width 8, height 12, dpi 600)## Figure 1c: Distal Uterine Tube Features for Cell Type Identification ##features - c(Dcn,Col1a1, # Fibroblasts Acta2,Myh11,Myl9, # Smooth MusclePdgfrb,Mcam,Cspg4, # PericyteSele,Vwf,Tek, # Blood EndothelialLyve1,Prox1,Icam1, # Lymph EndothelialEpcam,Krt8, # EpithelialFoxj1, # Ciliated EpithelialOvgp1, # Secretory EpithelialSlc1a3,Pax8,Cd44,Itga6, # Stem-like Epithelieal Ptprc, # Immune Cd8a,Cd4,Cd3e, # T-Cell Klrc1,Runx3, # T/NK CellKlrd1, # NK CellAif1,Cd68,Csf1r,Itgax, # MacrophageHbb-bs, Hba-a1, # ErythrocytesFras1,Rspo1,Lrrn4,Msln, # MesothelialCyp11a1,Bcat1,Fkbp5,Spp1,Prlr) # Luteal Cellsall_dp - DotPlot(object Distal_Named, # Seurat objectassay RNA, # Name of assay to use. Default is the active assayfeatures features, # List of features (select one from above or create a new one)# Colors to be used in the gradientcol.min 0, # Minimum scaled average expression threshold (everything smaller will be set to this)col.max 2.5, # Maximum scaled average expression threshold (everything larger will be set to this)dot.min 0, # The fraction of cells at which to draw the smallest dot (default is 0)dot.scale 9, # Scale the size of the pointsgroup.by NULL, # How the cells are going to be groupedsplit.by NULL, # Whether to split the data (if you fo this make sure you have a different color for each variable)scale TRUE, # Whether the data is scaledscale.by radius, # Scale the size of the points by size or radiusscale.min NA, # Set lower limit for scalingscale.max NA # Set upper limit for scaling
) labs(x NULL, y NULL)scale_color_viridis_c(optionF,begin.4,end0.9, direction -1)geom_point(aes(sizepct.exp), shape 21, colourblack, stroke0.6)#theme_linedraw()guides(x guide_axis(angle 90)) theme(axis.text.x element_text(size 14 , face italic))theme(axis.text.y element_text(size 14))scale_y_discrete(limits rev(levels(Distal_Named)))theme(legend.title element_text(size 14))ggsave(filename FIG1c_all_distal_dotplot.pdf, plot all_dp, width 18, height 10, dpi 600)图2 #### Figure 2: Characterization of distal epithelial cell states ######## Packages Load ####library(dplyr)
library(patchwork)
library(Seurat)
library(harmony)
library(ggplot2)
library(cowplot)
library(SoupX)
library(DoubletFinder)
library(data.table)
library(parallel)
library(tidyverse)
library(SoupX)
library(ggrepel)library(ggplot2)
library(gplots)
library(RColorBrewer)
library(viridisLite)
library(Polychrome)
library(circlize)
library(NatParksPalettes)library(monocle3)#### Load Distal Epithelial Dataset ####Epi_Filter - readRDS(file ../dataset/Distal_Epi_Cells.rds , refhook NULL)Epi_Named - RenameIdents(Epi_Filter, 0 Spdef Secretory, 1 Slc1a3 Stem/Progenitor, 2 Cebpdhigh/Foxj1- Progenitor,3 Ciliated 1, 4 Ciliated 2, 5 Pax8low/Prom1 Cilia-forming, 6 Fibroblast-like,7 Slc1a3med/Sox9 Cilia-forming,8 Selenop/Gstm2high Secretory)Epi_Namedactive.ident - factor(x Epi_Namedactive.ident, levels c( c(Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming, Fibroblast-like,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2)))EpiSecStemMarkers - FindMarkers(Epi_Named, ident.1 Spdef Secretory,ident.2 Slc1a3 Stem/Progenitor)
write.csv(EpiSecStemMarkers, file 20240319_Staining_Markers2.csv)#### Figure 2a: Epithelial Cells of the Distal Uterine Tube ####epi_umap - DimPlot(object Epi_Named, # Seurat object reduction umap, # Axes for the plot (UMAP, PCA, etc.) repel TRUE, # Whether to repel the cluster labelslabel FALSE, # Whether to have cluster labels cols c(#35EFEF, #1#00A1C6, #2#2188F7, #3#EA68E1, #4#59D1AF, #5#B20224, #6#F28D86, #7#A374B5, #8#9000C6), #9pt.size 0.6, # Size of each dot is (0.1 is the smallest)label.size 0.5) # Font size for labels # You can add any ggplot2 1customizations herelabs(title Colored by Cluster) # Plot titleNoLegend() labs(xUMAP_1,yUMAP_2)ggsave(filename Fig2a_epi_umap.pdf, plot epi_umap, width 15, height 12, dpi 600)#### Figure 2c: Distal Uterine Tube Features for Epithelial Cell State Identification ####distal_features - c(Krt8,Epcam,Slc1a3,Cd44,Sox9,Ovgp1,Sox17,Pax8, Egr1,Itga6, Bmpr1b,Rhoj, Klf6,Msln,Cebpd,Dpp6, Sec14l3, Fam161a,Prom1, Ly6a, Kctd8, Adam8,Dcn, Col1a1, Col1a2, Timp3, Pdgfra,Lgals1,Upk1a, Thrsp,Spdef,Lcn2,Selenop, Gstm2,Foxj1,Fam183b,Rgs22,Dnali1, Mt1 , Dynlrb2,Cdh1)epi_dp - DotPlot(object Epi_Named, # Seurat objectassay RNA, # Name of assay to use. Default is the active assayfeatures distal_features, # List of features (select one from above or create a new one)# Colors to be used in the gradientcol.min 0, # Minimum scaled average expression threshold (everything smaller will be set to this)col.max 2.5, # Maximum scaled average expression threshold (everything larger will be set to this)dot.min 0, # The fraction of cells at which to draw the smallest dot (default is 0)dot.scale 4, # Scale the size of the pointsgroup.by NULL, # How the cells are going to be groupedsplit.by NULL, # Whether to split the data (if you fo this make sure you have a different color for each variable)scale TRUE, # Whether the data is scaledscale.by radius, # Scale the size of the points by size or radiusscale.min NA, # Set lower limit for scalingscale.max NA ) # Set upper limit for scalinglabs(x NULL, # x-axis labely NULL)scale_color_viridis_c(optionF,begin.4,end0.9, direction -1)geom_point(aes(sizepct.exp), shape 21, colourblack, stroke0.6)#theme_linedraw()guides(x guide_axis(angle 90))theme(axis.text.x element_text(size 8 , face italic))theme(axis.text.y element_text(size 9))theme(legend.title element_text(size 9))theme(legend.text element_text(size 8)) scale_y_discrete(limits c(Ciliated 2,Ciliated 1,Selenop/Gstm2high Secretory,Spdef Secretory,Fibroblast-like,Pax8low/Prom1 Cilia-forming, Slc1a3med/Sox9 Cilia-forming,Cebpdhigh/Foxj1- Progenitor,Slc1a3 Stem/Progenitor))ggsave(filename Fig2b_epi_dot_plot.pdf, plot epi_dp, width 8.3, height 4.0625, dpi 600)#### Load Distal Epithelial Pseudotime Dataset ####Distal_PHATE - readRDS(file ../dataset/Distal_Epi_PHATE.rds , refhook NULL)Beeg_PHATE - Distal_PHATEBeeg_PHATEreductions[[phate]]cell.embeddings - Distal_PHATEreductions[[phate]]cell.embeddings*100cds - readRDS(file ../dataset/Distal_Epi_PHATE_Monocle3.rds , refhook NULL)#### Figure 2b: PHATE embedding for differentiation trajectory of distal epithelial cells ####phate_dif - DimPlot(Beeg_PHATE , reduction phate, cols c(#B20224, #35EFEF, #00A1C6, #A374B5, #9000C6, #EA68E1, #59D1AF, #2188F7, #F28D86),pt.size 0.7,shuffle TRUE,seed 0,label FALSE) labs(title Colored by Cluster) # Plot titleNoLegend() labs(xUMAP_1,yUMAP_2)ggsave(filename Fig3a_epi_phate.pdf, plot phate_dif, width 15, height 12, dpi 600)#### Figure 2d: PHATE and Monocle3 differentiation trajectory path ####pseudtotime - plot_cells(cds, color_cells_by pseudotime,label_cell_groupsFALSE,label_leavesFALSE,label_branch_pointsFALSE,graph_label_size0,cell_size .01,cell_stroke 1)theme(axis.title.x element_blank())theme(axis.title.y element_blank())theme(axis.line.x element_blank())theme(axis.line.y element_blank())theme(axis.ticks.x element_blank())theme(axis.ticks.y element_blank())theme(axis.text.x element_blank())theme(axis.text.y element_blank())theme(legend.text element_text(size 12))ggsave(filename Fig3b_epi_pseudtotime.pdf, plot pseudtotime, width 18, height 12, dpi 600)#### Figure 2e: Slc1a3 PHATE Feature Plot ####Slc1a3_PHATE - FeaturePlot(Beeg_PHATE, features c(Slc1a3), reduction phate, pt.size 0.5)scale_color_viridis_c(optionF,begin.4,end0.99, direction -1)theme(plot.title element_text(size 32,face bold.italic))theme(axis.title.x element_blank())theme(axis.title.y element_blank())theme(axis.line.x element_blank())theme(axis.line.y element_blank())theme(axis.ticks.x element_blank())theme(axis.ticks.y element_blank())theme(axis.text.x element_blank())theme(axis.text.y element_blank())theme(legend.text element_text(size 12))ggsave(filename Fig3c_SLC1A3_PHATE.pdf, plot Slc1a3_PHATE, width 18, height 9, dpi 600)#### Figure 2f: Pax8 PHATE Feature Plot ####Pax8_PHATE - FeaturePlot(Beeg_PHATE, features c(Pax8), reduction phate, pt.size 0.5)scale_color_viridis_c(optionF,begin.4,end0.99, direction -1)theme(plot.title element_text(size 32,face bold.italic))theme(axis.title.x element_blank())theme(axis.title.y element_blank())theme(axis.line.x element_blank())theme(axis.line.y element_blank())theme(axis.ticks.x element_blank())theme(axis.ticks.y element_blank())theme(axis.text.x element_blank())theme(axis.text.y element_blank())theme(legend.text element_text(size 12))ggsave(filename Fig3d_PAX8_PHATE.pdf, plot Pax8_PHATE, width 18, height 9, dpi 600)图6 #### Figure 6: Identification of cancer-prone cell states ######## Packages Load ####library(dplyr)
library(patchwork)
library(Seurat)
library(harmony)
library(ggplot2)
library(cowplot)
library(SoupX)
library(DoubletFinder)
library(data.table)
library(parallel)
library(tidyverse)
library(SoupX)
library(ggrepel)library(ggplot2)
library(gplots)
library(RColorBrewer)
library(viridisLite)
library(Polychrome)
library(circlize)
library(NatParksPalettes)library(monocle3)
library(ComplexHeatmap)
library(ggExtra)
library(gridExtra)
library(egg)library(scales)#### Load and Prepare Distal Epithelial and Epithelial Pseudotime Datasets ####Distal_PHATE - readRDS(file ../dataset/Distal_Epi_PHATE.rds , refhook NULL)cds - readRDS(file ../dataset/Distal_Epi_PHATE_Monocle3.rds , refhook NULL)Epi_Filter - readRDS(file ../dataset/Distal_Epi_Cells.rds , refhook NULL)Epi_Named - RenameIdents(Epi_Filter, 0 Spdef Secretory, 1 Slc1a3 Stem/Progenitor, 2 Cebpdhigh/Foxj1- Progenitor,3 Ciliated 1, 4 Ciliated 2, 5 Pax8low/Prom1 Cilia-forming, 6 Fibroblast-like,7 Slc1a3med/Sox9 Cilia-forming,8 Selenop/Gstm2high Secretory)Epi_Namedactive.ident - factor(x Epi_Namedactive.ident, levels c( c(Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming, Fibroblast-like,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2)))## Calculate Pseudotime Values ##pseudo - pseudotime(cds)Distal_PHATEmeta.data$Pseudotime - pseudo # Add to Seurat Metadata## Subset Seurat Object ##color_cells - DimPlot(Distal_PHATE , reduction phate, cols c(#B20224, #1#35EFEF, #2#00A1C6, #3#A374B5, #4#9000C6, #5#EA68E1, #6lightgrey, #7#2188F7, #8#F28D86),pt.size 0.7,shuffle TRUE,seed 0,label FALSE)## Psuedotime and Lineage Assignment ##cellID - rownames(Distal_PHATEreductions$phatecell.embeddings)
phate_embeddings - Distal_PHATEreductions$phatecell.embeddings
pseudotime_vals - Distal_PHATEmeta.data$Pseudotimecombined_data - data.frame(cellID, phate_embeddings, pseudotime_vals)# Calculate the Average PHATE_1 Value for Pseudotime Points 0 #
avg_phate_1 - mean(phate_embeddings[pseudotime_vals 0, 1])# Pseudotime Values lower than avge PHATE_1 Embedding will be Negative to split lineages
combined_data$Split_Pseudo - ifelse(phate_embeddings[, 1] avg_phate_1, -pseudotime_vals, pseudotime_vals)# Define Lineage #
combined_data$lineage - ifelse(combined_data$PHATE_1 avg_phate_1, Secretory,ifelse(combined_data$PHATE_1 avg_phate_1, Ciliogenic, Progenitor))Distal_PHATE$Pseudotime_Adj - combined_data$Split_Pseudo
Distal_PHATE$Lineage - combined_data$lineage# Subset #Pseudotime_Lineage - subset(Distal_PHATE, idents c(Secretory 1,Secretory 2,Msln Progenitor,Slc1a3/Sox9 Cilia-forming,Pax8/Prom1 Cilia-forming,Progenitor,Ciliated 1,Ciliated 2))## Set Bins ##bins - cut_number(Pseudotime_Lineagemeta.data$Pseudotime_Adj , 40) # Evenly distribute bins Pseudotime_Lineagemeta.data$Bin - bins # Metadata for Bins## Set Idents to PSeudoime Bin ##time_ident - SetIdent(Pseudotime_Lineage, value Pseudotime_Lineagemeta.data$Bin)av.exp - AverageExpression(time_ident, return.seurat T)$RNA # Calculate Avg log normalized expression# Calculates Average Expression for Each Bin #
# if you set return.seuratT, NormalizeData is called which by default performs log-normalization #
# Reported as avg log normalized expression ##### Figure 6c: PHATE embedding for differentiation trajectory of distal epithelial cells ##### Create the stacked barplot
rainbow20 - c(#FF0000,#FF6000,#FF8000,#FFA000,#FFC000,#FFE000,#FFFF00,#E0FF00,#C0FF00,#A0FF00,#00FF00,#00FFA0,#00F0FF,#00A0FF,#0020FF,#4000FF,#8000FF,#A000FF,#C000FF,#E000FF)rainbow_pseudo - DimPlot(Pseudotime_Lineage , reduction phate, cols c(rev(rainbow20),rainbow20),pt.size 1.2,shuffle TRUE,seed 0,label FALSE,group.by Bin) NoLegend()ggsave(filename rainbow_pseudo.pdf, plot rainbow_pseudo, width 20, height 10, dpi 600)#### Figure 6d: PHATE embedding for differentiation trajectory of distal epithelial cells ###### Pseudotime Scale Bar ##list - 1:40
colors c(rev(rainbow20),rainbow20)
df - data.frame(data list, color colors)pseudo_bar - ggplot(df, aes(x 1:40, y 1, fill color)) geom_bar(stat identity,position fill, color black, size 0, width 1) scale_fill_identity() theme_void() theme(axis.line element_blank(),axis.ticks element_blank(),axis.text element_blank(),axis.title element_blank())ggsave(filename pseudo_bar.pdf, plot pseudo_bar, width 0.98, height 0.19, dpi 600)## Plot gene list across pseudotime bin ##features - c(Upk1a, Spdef, Ovgp1, Gstm2, Selenop, Msln, Slc1a3,Itga6, Pax8,Ecrg4, Sox5, Pde4b, Lcn2,Klf6,Trp53 , Trp73,Krt5,Foxa2,Prom1,Clstn2,Spef2,Dnah12,Foxj1, Fam166c , Cfap126,Fam183b)# Create Bin List and expression of features #bin_list - unique(Pseudotime_Lineagemeta.data$Bin) plot_info - as.data.frame(av.exp[features, ]) # Call Avg Expression for featuresz_score - transform(plot_info, SDapply(plot_info,1, mean, na.rm TRUE))
z_score - transform(z_score, MEANapply(plot_info,1, sd, na.rm TRUE))z_score1 - (plot_info-z_score$MEAN)/z_score$SDplot_info$y - rownames(plot_info) # y values as features
z_score1$y - rownames(plot_info)plot_info - gather(data plot_info, x, expression, bin_list) #set plot
z_score1 - gather(data z_score1, x, z_score, bin_list) #set plot# Create Cell Clusters DF #Labeled_Pseudotime_Lineage - RenameIdents(Pseudotime_Lineage, Secretory 1 Spdef Secretory, Progenitor Slc1a3 Stem/Progenitor, Msln Progenitor Cebpdhigh/Foxj1- Progenitor,Ciliated 1 Ciliated 1, Ciliated 2 Ciliated 2, Pax8/Prom1 Cilia-forming Pax8low/Prom1 Cilia-forming, Fibroblast-like Fibroblast-like, #removedSlc1a3/Sox9 Cilia-forming Slc1a3med/Sox9 Cilia-forming,Secretory 2 Selenop/Gstm2high Secretory)cluster_table - table(Labeled_Pseudotime_Lineageactive.ident, Labeled_Pseudotime_Lineagemeta.data$Bin)clusters - data.frame(cluster_table)clusters - clusters %% group_by(Var2) %%mutate(Perc Freq / sum(Freq))# Create Pseudotime DF #pseudotime_table - table(seq(1, length(bin_list), 1), unique(Labeled_Pseudotime_Lineagemeta.data$Bin),seq(1, length(bin_list), 1))pseudotime_bins - data.frame(pseudotime_table) # calculate max and min z-scores
max_z - max(z_score1$z_score, na.rm TRUE)
min_z - min(z_score1$z_score, na.rm TRUE)# set color for outliers
outlier_color - ifelse(z_score1$z_score max_z | z_score1$z_score min_z, ifelse(z_score1$z_score 0, #AD1F24, #51A6DC), #e2e2e2)## Plot Gene Expression ### Set different na.value options for positive and negative values
na_color_pos - #AD1F24 # color for positive NA values
na_color_neg - #51A6DC # color for negative NA valuescustom_bin_names - c(paste0(S, 20:1), paste0(C, 1:20))figure - ggplot(z_score1, aes(x, y, fill z_score)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradientn(colorsc(#1984c5, #e2e2e2, #c23728), name Average Expression \nZ-Score, limits c(-3,3), na.value ifelse(is.na(z_score1) z_score1 0, na_color_pos, ifelse(is.na(z_score1) z_score1 0, na_color_neg, grey50)),oob scales::squish)scale_x_discrete(limits sort(bin_list) , labels custom_bin_names)scale_y_discrete(limits rev(features))theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold), # Text size throughout the plotaxis.text.x element_text(color black, angle 0, hjust 0.5, size 10, face bold), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold.italic))theme(plot.title element_blank(),plot.marginunit(c(-0.5,1,1,1), cm))## Plot Cluster Percentage ##Spdef Secretory - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Spdef Secretory)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(1,1,1,1), cm))Selenop/Gstm2high Secretory - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Selenop/Gstm2high Secretory)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))Cebpdhigh/Foxj1- Progenitor - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Cebpdhigh/Foxj1- Progenitor)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))Slc1a3 Stem/Progenitor - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Slc1a3 Stem/Progenitor)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))Slc1a3med/Sox9 Cilia-forming - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Slc1a3med/Sox9 Cilia-forming)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))Pax8low/Prom1 Cilia-forming - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Pax8low/Prom1 Cilia-forming)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))Ciliated 1 - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Ciliated 1)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))Ciliated 2 - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Ciliated 2)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))## Plot Pseudotime Color ##list - 1:40
colors c(rev(rainbow20),rainbow20)
df - data.frame(data list, color colors)binning - ggplot(df, aes(x 1:40, y 1, fill color)) geom_bar(stat identity,position fill, color black, size 1, width 1) scale_fill_identity() theme_void() theme(axis.line element_blank(),axis.ticks element_blank(),axis.text element_blank(),axis.title element_blank())scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Pseudotime Bin )theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, vjust .75, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))### Combine Plots ###psuedotime_lineage - ggarrange(Spdef Secretory,Selenop/Gstm2high Secretory,Cebpdhigh/Foxj1- Progenitor,Slc1a3 Stem/Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming,Ciliated 1,Ciliated 2,binning,figure , ncol1,heights c(2, 2, 2, 2, 2, 2, 2, 2, 2, (2*length(features)),widths c(3)),padding unit(0.01))ggsave(filename FIG6d_psuedotime_lineage.pdf, plot psuedotime_lineage, width 18, height 9, dpi 600)#### Figure 6e: Stacked violin plots for cancer-prone cell states ####### Stacked Violin Plot Function ####https://divingintogeneticsandgenomics.rbind.io/post/stacked-violin-plot-for-visualizing-single-cell-data-in-seurat/## remove the x-axis text and tick
## plot.margin to adjust the white space between each plot.
## ... pass any arguments to VlnPlot in Seurat
modify_vlnplot - function(obj, feature, pt.size 0, plot.margin unit(c(-0.75, 0, -0.75, 0), cm),...) {p- VlnPlot(obj, features feature, pt.size pt.size, ... ) xlab() ylab(feature) ggtitle() theme(legend.position none, axis.text.x element_blank(), axis.ticks.x element_blank(), axis.title.y element_text(size rel(1), angle 0, face bold.italic), axis.text.y element_text(size rel(1)), plot.margin plot.margin ) return(p)
}## extract the max value of the y axis
extract_max- function(p){ymax- max(ggplot_build(p)$layout$panel_scales_y[[1]]$range$range)return(ceiling(ymax))
}## main function
StackedVlnPlot- function(obj, features,pt.size 0, plot.margin unit(c(-0.75, 0, -0.75, 0), cm),...) {plot_list- purrr::map(features, function(x) modify_vlnplot(obj obj,feature x, ...))# Add back x-axis title to bottom plot. patchwork is going to support this?plot_list[[length(plot_list)]]- plot_list[[length(plot_list)]] theme(axis.text.xelement_text(angle 60, hjust1, vjust0.95), axis.ticks.x element_line())# change the y-axis tick to only max value ymaxs- purrr::map_dbl(plot_list, extract_max)plot_list- purrr::map2(plot_list, ymaxs, function(x,y) x scale_y_continuous(breaks c(y)) expand_limits(y y))p- patchwork::wrap_plots(plotlist plot_list, ncol 1)return(p)
}features- c(Slc1a3, Pax8 , Trp73 , Prom1 )features- c(Pax8, Ovgp1 , Lcn2 , Upk1a , Spdef ,Thrsp )colors - c(#35EFEF, #1#00A1C6, #2#2188F7, #3#EA68E1, #4##59D1AF, #5#B20224, #6#F28D86, #7#A374B5, #8#9000C6)
No_Fibro - subset(x Epi_Named, idents c(Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming, #Fibroblast-like,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2))stack_vln - StackedVlnPlot(obj No_Fibro, features features, slot data,pt.size 0,cols c(#35EFEF, #1#00A1C6, #2#2188F7, #3#EA68E1, #4##59D1AF, #5#B20224, #6#F28D86, #7#A374B5, #8#9000C6)) #9theme(plot.title element_text(size 32, face bold.italic))scale_x_discrete(limits c(Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming, #Fibroblast-like,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2))theme(axis.text.x element_text(size 16, angle 60))theme(axis.text.y element_text(size 14))theme(axis.title.y.left element_text(size 16))ggsave(filename FIG6e_stacked_vln_noFibro.pdf, plot stack_vln, width 18, height 12, dpi 600)#### Figure 6f: Krt5 expression within epithelial cell states ####ggsave(filename Krt5_dp_others.pdf, plot Krt5_dp_others, width 1.89*8, height 3.06*8, dpi 600)Krt5_dp - DotPlot(object No_Fibro, # Seurat objectassay RNA, # Name of assay to use. Default is the active assayfeatures Krt5, # List of features (select one from above or create a new one)# Colors to be used in the gradientcol.min 0, # Minimum scaled average expression threshold (everything smaller will be set to this)col.max 2.5, # Maximum scaled average expression threshold (everything larger will be set to this)dot.min 0, # The fraction of cells at which to draw the smallest dot (default is 0)dot.scale 24, # Scale the size of the pointsgroup.by NULL, # How the cells are going to be groupedsplit.by NULL, # Whether to split the data (if you fo this make sure you have a different color for each variable)scale TRUE, # Whether the data is scaledscale.by radius, # Scale the size of the points by size or radiusscale.min NA, # Set lower limit for scalingscale.max NA ) # Set upper limit for scalinglabs(x NULL, # x-axis labely NULL)scale_color_viridis_c(optionF,begin.4,end0.9, direction -1)geom_point(aes(sizepct.exp), shape 21, colourblack, stroke0.7)theme_linedraw(base_line_size 5)guides(x guide_axis(angle 90))theme(axis.text.x element_text(size 32 , face italic))theme(axis.text.y element_text(size 32))theme(legend.title element_text(size 12)) scale_y_discrete(limits c(Ciliated 2,Ciliated 1,Selenop/Gstm2high Secretory,Spdef Secretory,#Fibroblast-like,Pax8low/Prom1 Cilia-forming, Slc1a3med/Sox9 Cilia-forming,Cebpdhigh/Foxj1- Progenitor,Slc1a3 Stem/Progenitor))ggsave(filename Krt5_dp_noFibro.pdf, plot Krt5_dp, width 1.89*8, height 3.06*8, dpi 600)附图2 #### Figure Supp 2: Characterization of distal epithelial cell states ######## Packages Load ####library(dplyr)
library(patchwork)
library(Seurat)
library(harmony)
library(ggplot2)
library(cowplot)
library(SoupX)
library(DoubletFinder)
library(data.table)
library(parallel)
library(tidyverse)
library(SoupX)
library(ggrepel)library(ggplot2)
library(gplots)
library(RColorBrewer)
library(viridisLite)
library(Polychrome)
library(circlize)
library(NatParksPalettes)#### Unprocessed and Processed Distal Dataset ####All.merged - readRDS(file ../dataset/Unfiltered_Mouse_Distal.rds, refhook NULL) # Prior to Quality ControlDistal - readRDS(file ../dataset/Distal_Filtered_Cells.rds , refhook NULL) # After to Quality ControlDistal_Named - RenameIdents(Distal, 0 Fibroblast 1, 1 Fibroblast 2, 2 Secretory Epithelial,3 Smooth Muscle, 4 Ciliated Epithelial 1, 5 Fibroblast 3, 6 Stem-like Epithelial 1,7 Stem-like Epithelial 2,8 Ciliated Epithelial 2, 9 Blood Endothelial, 10 Pericyte, 11 Intermediate Epithelial, 12 T/NK Cell, 13 Epithelial/Fibroblast, 14 Macrophage, 15 Erythrocyte, 16 Luteal,17 Mesothelial,18 Lymph Endothelial/Epithelial) # Remove cluster due few data points and suspected doubletDistal_Namedactive.ident - factor(x Distal_Namedactive.ident, levels c(Fibroblast 1,Fibroblast 2,Fibroblast 3,Smooth Muscle,Pericyte,Blood Endothelial,Lymph Endothelial/Epithelial,Epithelial/Fibroblast,Stem-like Epithelial 1,Stem-like Epithelial 2,Intermediate Epithelial,Secretory Epithelial,Ciliated Epithelial 1,Ciliated Epithelial 2,T/NK Cell,Macrophage,Erythrocyte,Mesothelial,Luteal))Distal_Named - subset(Distal_Named, idents c(Fibroblast 1,Fibroblast 2,Fibroblast 3,Smooth Muscle,Pericyte,Blood Endothelial,Epithelial/Fibroblast,Stem-like Epithelial 1,Stem-like Epithelial 2,Intermediate Epithelial,Secretory Epithelial,Ciliated Epithelial 1,Ciliated Epithelial 2,T/NK Cell,Macrophage,Erythrocyte,Mesothelial,Luteal))Distal_Named - SetIdent(Distal_Named, value Distal_Namedactive.ident)Epi_Filter - readRDS(file ../dataset/Distal_Epi_Cells.rds , refhook NULL)Epi_Named - RenameIdents(Epi_Filter, 0 Spdef Secretory, 1 Slc1a3 Stem/Progenitor, 2 Cebpdhigh/Foxj1- Progenitor,3 Ciliated 1, 4 Ciliated 2, 5 Pax8low/Prom1 Cilia-forming, 6 Fibroblast-like,7 Slc1a3med/Sox9 Cilia-forming,8 Selenop/Gstm2high Secretory)Epi_Namedactive.ident - factor(x Epi_Namedactive.ident, levels c( c(Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming, Fibroblast-like,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2)))#### Figure Supp 2a: Unfilitered % MT Genes ####unfiltered_MT - VlnPlot(All.merged, features c(percent.mt), group.by Sample, pt.size 0,cols natparks.pals(nameArches2,n3))theme(legend.position none)theme(axis.text.x element_text(size 16)) # Change X-Axis Text Sizetheme(axis.text.y element_text(size 16)) # Change Y-Axis Text Sizetheme(axis.title.y element_text(size 18)) # Change Y-Axis Title Text Sizetheme(plot.title element_text(size 32,face bold.italic))theme(axis.title.x element_blank())ggsave(filename FIGs2a_unfiltered_MT.pdf, plot unfiltered_MT, width 12, height 9, dpi 600)#exprs - as.data.frame(FetchData(object All.merged, vars c(nCount_RNA , Sample)))#df_new - filter(exprs, Sample mD1)#mean(df_new$nCount_RNA)
#sd(df_new$nCount_RNA)# Remove the original Value column if needed
df_new - df_new %%select(-Value)x - spread(exprs, Sample, percent.mt )
mean(exprs)
#### Figure Supp 2b: Unfilitered nFeature RNA #### unfiltered_nFeature - VlnPlot(All.merged, features c(nFeature_RNA), group.by Sample, pt.size 0,cols natparks.pals(nameArches2,n3))theme(legend.position none)theme(axis.text.x element_text(size 16))theme(axis.text.y element_text(size 16))theme(axis.title.y element_text(size 18))theme(plot.title element_text(size 32,face bold.italic))theme(axis.title.x element_blank()) # Change object to visualize other samplesggsave(filename FIGs2b_unfiltered_nFeature.pdf, plot unfiltered_nFeature, width 12, height 9, dpi 600)#### Figure Supp 2c: Unfilitered nCount RNA #### unfiltered_nCount - VlnPlot(All.merged, features c(nCount_RNA), group.by Sample, pt.size 0,cols natparks.pals(nameArches2,n3))theme(legend.position none)theme(axis.text.x element_text(size 16))theme(axis.text.y element_text(size 16))theme(axis.title.y element_text(size 18))theme(plot.title element_text(size 32,face bold.italic))theme(axis.title.x element_blank()) # Change object to visualize other samplesggsave(filename FIGs2c_unfiltered_nCount.pdf, plot unfiltered_nCount, width 12, height 9, dpi 600)#### Figure Supp 2d: Doublets All Cells #### All_Doublet - DimPlot(object Distal, reduction umap, group.by Doublet,cols c( #ffb6c1, #380b11),repel TRUE, label F, pt.size 1.2, order c(Doublet,Singlet),label.size 5) labs(xUMAP_1,yUMAP_2)ggsave(filename FIGs2d1_All_Doublet_umap.pdf, plot All_Doublet, width 22, height 17, dpi 600)## Stacked Bar Doublets ##table - table(Distal_Namedactive.ident ,Distal_Namedmeta.data$Doublet) # Create a table of countsdf - data.frame(table) doublet - ggplot(data df, # Dataset to use for plot. Needs to be a data.frame aes(x Var1, # Variable to plot on the x-axisy Freq, # Variable to plot on the y-axisfill factor(Var2, # Variable to fill the barslevels c(Doublet,Singlet)))) # Order of the stacked barstheme_classic() # ggplot2 theme# Bar plotgeom_bar(position fill, # Position of bars. Fill means the bars are stacked.stat identity, # Height of bars represent values in the datasize 1) # Size of bars# Color schemescale_fill_manual(Doublet, limits c(Doublet,Singlet),values c(#8B0000,#808080)) # Add plot labelslabs(x NULL, # x-axis labely Fraction of Cells) # y-axis labeltheme(text element_text(size 15), # Text size throughout the plotaxis.text.x element_text(color black, angle 60, hjust 1, size 11), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1)) # Text color and horizontal adjustment on y-axisscale_x_discrete(limits (c(Intermediate Epithelial,Epithelial/Fibroblast,Stem-like Epithelial 1,Ciliated Epithelial 1,Erythrocyte,Smooth Muscle,Stem-like Epithelial 2,Mesothelial,Blood Endothelial,Pericyte,Fibroblast 2,Secretory Epithelial,Fibroblast 1,Ciliated Epithelial 2,Fibroblast 3,T/NK Cell,Macrophage,Luteal)))coord_flip()ggsave(filename FIGs2d2_doublet_quant.pdf, plot doublet, width 10, height 16, dpi 600)#### Figure Supp 2e: Sample Distribution for All Cells #### table - table(Distal_Namedactive.ident ,Distal_Namedmeta.data$Sample) # Create a table of countsdf - data.frame(table) table2 - table(Epi_Namedmeta.data$Sample)all_sample_dist - ggplot(data df, # Dataset to use for plot. Needs to be a data.frame aes(x Var1, # Variable to plot on the x-axisy Freq, # Variable to plot on the y-axisfill factor(Var2, # Variable to fill the barslevels c(mD1,mD2,mD4)))) # Order of the stacked barstheme_classic() # ggplot2 theme# Bar plotgeom_bar(position fill, # Position of bars. Fill means the bars are stacked.stat identity, # Height of bars represent values in the datasize 1) # Size of bars# Color schemescale_fill_manual(Location, limits c(mD1,mD2,mD4),values c(natparks.pals(nameArches2,n3))) # Add plot labelslabs(x NULL, # x-axis labely Fraction of Cells) # y-axis labeltheme(text element_text(size 15), # Text size throughout the plotaxis.text.x element_text(color black, angle 60, hjust 1, size 11), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1)) # Text color and horizontal adjustment on y-axisggsave(filename FIGs2f_all_sample_dist.pdf, plot epi_sample_dist, width 16, height 12, dpi 600)#### Figure Supp 2f: Doublets Epithelial Cells #### Epi_Doublet - DimPlot(object Epi_Named, reduction umap, group.by Doublet,cols c( #ffb6c1, #380b11),repel TRUE, label F, pt.size 1.2, order c(Doublet,Singlet),label.size 5) labs(xUMAP_1,yUMAP_2)ggsave(filename FIGs2e1_Epi_Doublet_umap.pdf, plot All_Doublet, width 22, height 17, dpi 600)## Stacked Bar Doublets ##table - table(Epi_Namedactive.ident ,Epi_Namedmeta.data$Doublet) # Create a table of countsdf - data.frame(table) epi_doublet - ggplot(data df, # Dataset to use for plot. Needs to be a data.frame aes(x Var1, # Variable to plot on the x-axisy Freq, # Variable to plot on the y-axisfill factor(Var2, # Variable to fill the barslevels c(Doublet,Singlet)))) # Order of the stacked barstheme_classic() # ggplot2 theme# Bar plotgeom_bar(position fill, # Position of bars. Fill means the bars are stacked.stat identity, # Height of bars represent values in the datasize 1) # Size of bars# Color schemescale_fill_manual(Doublet, limits c(Doublet,Singlet),values c(#8B0000,#808080)) # Add plot labelslabs(x NULL, # x-axis labely Fraction of Cells) # y-axis labeltheme(text element_text(size 15), # Text size throughout the plotaxis.text.x element_text(color black, angle 60, hjust 1, size 11), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1)) # Text color and horizontal adjustment on y-axisscale_x_discrete(limits (c(Slc1a3med/Sox9 Cilia-forming,Fibroblast-like,Pax8low/Prom1 Cilia-forming,Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Ciliated 1,Ciliated 2, Spdef Secretory,Selenop/Gstm2high Secretory)))coord_flip()ggsave(filename FIGs2e2_epi_doublet_quant.pdf, plot epi_doublet, width 10, height 16, dpi 600)#### Figure Supp 2g: Sample Distribution for Epithelial Cells #### table - table(Epi_Namedactive.ident ,Epi_Namedmeta.data$Sample) # Create a table of countsdf - data.frame(table) table2 - table(Epi_Namedmeta.data$Samlpe)epi_sample_dist - ggplot(data df, # Dataset to use for plot. Needs to be a data.frame aes(x Var1, # Variable to plot on the x-axisy Freq, # Variable to plot on the y-axisfill factor(Var2, # Variable to fill the barslevels c(mD1,mD2,mD4)))) # Order of the stacked barstheme_classic() # ggplot2 theme# Bar plotgeom_bar(position fill, # Position of bars. Fill means the bars are stacked.stat identity, # Height of bars represent values in the datasize 1) # Size of bars# Color schemescale_fill_manual(Location, limits c(mD1,mD2,mD4),values c(natparks.pals(nameArches2,n3))) # Add plot labelslabs(x NULL, # x-axis labely Fraction of Cells) # y-axis labeltheme(text element_text(size 15), # Text size throughout the plotaxis.text.x element_text(color black, angle 60, hjust 1, size 11), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1)) # Text color and horizontal adjustment on y-axisggsave(filename FIGs2g_epi_sample_dist.pdf, plot epi_sample_dist, width 16, height 12, dpi 600)#### Figure Supp 2h: Distal Tile Mosaic ####library(treemap)dist_cell_types - table(Idents(Distal_Named), Distal_Named$orig.ident)
dist_cell_type_df - as.data.frame(dist_cell_types)## Colors ##Fibroblasts - c(#FF9D00 , #FFB653 , #FFCB9A) # Oranges
FiboEpi - #FFE0B3 # Reddish Brown
Muscle - c(#E55451 , #FFB7B2) # Reds
Endothelial - c(#A0E6FF) # Reds
Epi -c(#6E3E6E,#8A2BE2,#604791,#CCCCFF,#DA70D6,#DF73FF) # Blues/Purples
Immune - c( #5A5E6B , #B8C2CC , #FC86AA) # Yellowish Brown
Meso - #9EFFFF # Pink
Lut - #9DCC00 # Greencolors - c(Fibroblasts, FiboEpi, Muscle, Endothelial, Epi, Immune, Meso, Lut)## Tile Mosaic ##distal_treemap - treemap(dist_cell_type_df, index Var1, vSize Freq, vColor colors, palette colors)ggsave(filename 20240612_all_distal_tile.pdf, plot distal_treemap, width 12, height 8, dpi 600)#### Figure Supp 2i: Epi Markers All Distal Cells ####### Stacked Violin Plot Function ####https://divingintogeneticsandgenomics.rbind.io/post/stacked-violin-plot-for-visualizing-single-cell-data-in-seurat/## remove the x-axis text and tick
## plot.margin to adjust the white space between each plot.
## ... pass any arguments to VlnPlot in Seuratmodify_vlnplot - function(obj, feature, pt.size 0, plot.margin unit(c(-0.75, 0, -0.75, 0), cm),...) {p- VlnPlot(obj, features feature, pt.size pt.size, ... ) xlab() ylab(feature) ggtitle() theme(legend.position none, axis.text.x element_blank(), axis.ticks.x element_blank(), axis.title.y element_text(size rel(1), angle 0, face bold.italic), axis.text.y element_text(size rel(1)), plot.margin plot.margin ) return(p)
}## extract the max value of the y axis
extract_max- function(p){ymax- max(ggplot_build(p)$layout$panel_scales_y[[1]]$range$range)return(ceiling(ymax))
}## main function
StackedVlnPlot- function(obj, features,pt.size 0, plot.margin unit(c(-0.75, 0, -0.75, 0), cm),...) {plot_list- purrr::map(features, function(x) modify_vlnplot(obj obj,feature x, ...))# Add back x-axis title to bottom plot. patchwork is going to support this?plot_list[[length(plot_list)]]- plot_list[[length(plot_list)]] theme(axis.text.xelement_text(angle 60, hjust1, vjust0.95), axis.ticks.x element_line())# change the y-axis tick to only max value ymaxs- purrr::map_dbl(plot_list, extract_max)# plot_list- purrr::map2(plot_list, ymaxs, function(x,y) x plot_list- purrr::map2(plot_list, c(5,5,8,5), function(x,y) x scale_y_continuous(breaks c(y)) expand_limits(y y))p- patchwork::wrap_plots(plotlist plot_list, ncol 1)return(p)
}features- c(Epcam, Krt8 , Ovgp1 , Foxj1 )Fibroblasts - c(#FF9D00 , #FFB653 , #FFCB9A) # Oranges
Muscle - c(#E55451 , #FFB7B2) # Reds
Endothelial - c(#A0E6FF) # Reds
FiboEpi - #FFE0B3 # Reddish Brown
Epi -c(#6E3E6E,#8A2BE2,#604791,#CCCCFF,#DA70D6,#DF73FF) # Blues/Purples
Immune - c( #5A5E6B , #B8C2CC , #FC86AA) # Yellowish Brown
Meso - #9EFFFF # Pink
Lut - #9DCC00 # Greencolors - c(Fibroblasts, FiboEpi, Muscle, Endothelial, Epi, Immune, Meso, Lut)stack_vln - StackedVlnPlot(obj Distal_Named, features features, slot data,pt.size 0,cols colors) #9theme(plot.title element_text(size 32, face bold.italic))scale_x_discrete(limits c(Fibroblast 1,Fibroblast 2,Fibroblast 3,Epithelial/Fibroblast,Smooth Muscle,Pericyte,Blood Endothelial,Stem-like Epithelial 1,Stem-like Epithelial 2,Intermediate Epithelial,Secretory Epithelial,Ciliated Epithelial 1,Ciliated Epithelial 2,T/NK Cell,Macrophage,Erythrocyte,Mesothelial,Luteal))theme(axis.text.x element_text(size 16, angle 60))theme(axis.text.y element_text(size 14))theme(axis.title.y.left element_text(size 16))ggsave(filename 20240612_All_Distal_stacked_vln.pdf, plot stack_vln, width 18, height 12, dpi 600)#### Figure Supp 2j: Epi Markers Distal Epi Cells ####Epi_Sub - subset(Epi_Named, idents c(Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming, Fibroblast-like,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2))colors - c(#35EFEF, #1#00A1C6, #2#2188F7, #3#EA68E1, #4#59D1AF, #5#B20224, #6#F28D86, #7#A374B5, #8#9000C6)stack_vln - StackedVlnPlot(obj Epi_Named, features features, slot data,pt.size 0,cols c(#35EFEF, #1#00A1C6, #2#2188F7, #3#EA68E1, #4#59D1AF, #5#B20224, #6#F28D86, #7#A374B5, #8#9000C6)) #9theme(plot.title element_text(size 32, face bold.italic))scale_x_discrete(limits c(Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming, Fibroblast-like,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2))theme(axis.text.x element_text(size 16, angle 60))theme(axis.text.y element_text(size 14))theme(axis.title.y.left element_text(size 16))ggsave(filename 20240612_Distal_Epi_stacked_vln.pdf, plot stack_vln, width 18, height 12, dpi 600)附图3
#### Figure Supp 3: Doublet detection of fibroblast and epithelial markers ######## Packages Load ####library(dplyr)
library(patchwork)
library(Seurat)
library(harmony)
library(ggplot2)
library(cowplot)
library(SoupX)
library(DoubletFinder)
library(data.table)
library(parallel)
library(tidyverse)
library(SoupX)
library(ggrepel)library(ggplot2)
library(gplots)
library(RColorBrewer)
library(viridisLite)
library(Polychrome)
library(circlize)
library(NatParksPalettes)#### Load Distal Epithelial Dataset ####Epi_Filter - readRDS(file ../dataset/Distal_Epi_Cells.rds , refhook NULL)Epi_Named - RenameIdents(Epi_Filter, 0 Spdef Secretory, 1 Slc1a3 Stem/Progenitor, 2 Cebpdhigh/Foxj1- Progenitor,3 Ciliated 1, 4 Ciliated 2, 5 Pax8low/Prom1 Cilia-forming, 6 Fibroblast-like,7 Slc1a3med/Sox9 Cilia-forming,8 Selenop/Gstm2high Secretory)Epi_Namedactive.ident - factor(x Epi_Namedactive.ident, levels c( c(Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming, Fibroblast-like,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2)))table(Epi_Namedactive.ident)#### Plot Compilation ####feature_scatter - FeatureScatter( Fibroblast, Krt8,Col1a1,cols c(#35EFEF, #1#00A1C6, #2#2188F7, #3#EA68E1, #4#59D1AF, #5#B20224, #6#F28D86, #7#A374B5, #8#9000C6))x - DotPlot(Epi_Named , features c(Krt8 , Col1a1))
write.csv(x$data , doublet_data.csv)Fibroblast - subset(Epi_Named, idents c(Fibroblast-like))custom_labels - function(x) {ifelse(x %% 1 0, as.character(x), )
}feature_scatter - FeatureScatter( Fibroblast, Krt8,Col1a1,cols black) # Scatter plotNoLegend()labs(title NULL)theme(panel.grid.major element_line(color grey, size 0.5),panel.grid.minor element_blank())theme(axis.text.x element_text(color black, size 12),axis.title.x element_text(color black, size 14, face bold.italic),axis.ticks.x element_line(color black),axis.text.y element_text(color black, size 12),axis.title.y element_text(color black, size 14, face bold.italic),axis.ticks.y element_line(color black),plot.margin unit(c(0, 0, 0, 0), mm))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on x-axisscale_y_continuous(limits c(0,5) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on y-axisplot_data - as.data.frame(feature_scatter$data)# Create density plot for x-axis
density_x - ggplot(plot_data, aes(x Krt8 , fill black)) geom_density() theme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_blank(),axis.title.x element_blank(),axis.ticks.x element_blank(),axis.text.y element_text(color white, size 12),axis.title.y element_text(color white, size 14, face bold.italic),axis.ticks.y element_line(color white),plot.margin unit(c(0, 0, 0, 0), mm))scale_y_continuous(labels function(y) sprintf(%.0f, y)) NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on x-axis# Create density plot for y-axis
density_y - ggplot(plot_data, aes(x Col1a1 , fill black)) geom_density() coord_flip() # Flip axes for y-density plottheme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_text(color white, size 12),axis.title.x element_text(color white, size 14, face bold.italic),axis.ticks.x element_line(color white),axis.text.y element_blank(),axis.title.y element_blank(),axis.ticks.y element_blank(),plot.margin unit(c(0, 0, 0, 0), mm))NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,5))scale_y_continuous(limits c(0,1) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on y-axis# Arrange plotstop_row - cowplot::plot_grid(density_x, NULL,nrow 1, rel_widths c(3, 1), rel_heights c(0.5,0.5))bottom_row - cowplot::plot_grid(feature_scatter,density_y,nrow 1, rel_widths c(3, 1), rel_heights c(3,3))combined_plot - cowplot::plot_grid(top_row , bottom_row,nrow 2 , rel_widths c(3, 1), rel_heights c(0.25,3))ggsave(filename 20240221_Fibroblast_Doublet.pdf, plot combined_plot, width 12, height 12, dpi 600)## Stem ##c(Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2)Stem - subset(Epi_Named, idents c(Slc1a3 Stem/Progenitor))custom_labels - function(x) {ifelse(x %% 1 0, as.character(x), )
}feature_scatter - FeatureScatter( Stem, Krt8,Col1a1,cols black) # Scatter plotNoLegend()labs(title NULL)theme(panel.grid.major element_line(color grey, size 0.5),panel.grid.minor element_blank())theme(axis.text.x element_text(color black, size 12),axis.title.x element_text(color black, size 14, face bold.italic),axis.ticks.x element_line(color black),axis.text.y element_text(color black, size 12),axis.title.y element_text(color black, size 14, face bold.italic),axis.ticks.y element_line(color black),plot.margin unit(c(0, 0, 0, 0), mm))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on x-axisscale_y_continuous(limits c(0,5) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on y-axisplot_data - as.data.frame(feature_scatter$data)# Create density plot for x-axis
density_x - ggplot(plot_data, aes(x Krt8 , fill black)) geom_density() theme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_blank(),axis.title.x element_blank(),axis.ticks.x element_blank(),axis.text.y element_text(color white, size 12),axis.title.y element_text(color white, size 14, face bold.italic),axis.ticks.y element_line(color white),plot.margin unit(c(0, 0, 0, 0), mm))scale_y_continuous(labels function(y) sprintf(%.0f, y)) NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on x-axis# Create density plot for y-axis
density_y - ggplot(plot_data, aes(x Col1a1 , fill black)) geom_density() coord_flip() # Flip axes for y-density plottheme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_text(color white, size 12),axis.title.x element_text(color white, size 14, face bold.italic),axis.ticks.x element_line(color white),axis.text.y element_blank(),axis.title.y element_blank(),axis.ticks.y element_blank(),plot.margin unit(c(0, 0, 0, 0), mm))NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,5))scale_y_continuous(limits c(0,1) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on y-axis# Arrange plotstop_row - cowplot::plot_grid(density_x, NULL,nrow 1, rel_widths c(3, 1), rel_heights c(0.5,0.5))bottom_row - cowplot::plot_grid(feature_scatter,density_y,nrow 1, rel_widths c(3, 1), rel_heights c(3,3))combined_plot - cowplot::plot_grid(top_row , bottom_row,nrow 2 , rel_widths c(3, 1), rel_heights c(0.25,3))ggsave(filename 20240221_Slc1a3_Doublet.pdf, plot combined_plot, width 12, height 12, dpi 600)## Prog ##c(Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2)Prog - subset(Epi_Named, idents c(Cebpdhigh/Foxj1- Progenitor))custom_labels - function(x) {ifelse(x %% 1 0, as.character(x), )
}feature_scatter - FeatureScatter( Prog, Krt8,Col1a1,cols black) # Scatter plotNoLegend()labs(title NULL)theme(panel.grid.major element_line(color grey, size 0.5),panel.grid.minor element_blank())theme(axis.text.x element_text(color black, size 12),axis.title.x element_text(color black, size 14, face bold.italic),axis.ticks.x element_line(color black),axis.text.y element_text(color black, size 12),axis.title.y element_text(color black, size 14, face bold.italic),axis.ticks.y element_line(color black),plot.margin unit(c(0, 0, 0, 0), mm))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on x-axisscale_y_continuous(limits c(0,5) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on y-axisplot_data - as.data.frame(feature_scatter$data)# Create density plot for x-axis
density_x - ggplot(plot_data, aes(x Krt8 , fill black)) geom_density() theme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_blank(),axis.title.x element_blank(),axis.ticks.x element_blank(),axis.text.y element_text(color white, size 12),axis.title.y element_text(color white, size 14, face bold.italic),axis.ticks.y element_line(color white),plot.margin unit(c(0, 0, 0, 0), mm))scale_y_continuous(labels function(y) sprintf(%.0f, y)) NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on x-axis# Create density plot for y-axis
density_y - ggplot(plot_data, aes(x Col1a1 , fill black)) geom_density() coord_flip() # Flip axes for y-density plottheme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_text(color white, size 12),axis.title.x element_text(color white, size 14, face bold.italic),axis.ticks.x element_line(color white),axis.text.y element_blank(),axis.title.y element_blank(),axis.ticks.y element_blank(),plot.margin unit(c(0, 0, 0, 0), mm))NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,5))scale_y_continuous(limits c(0,1) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on y-axis# Arrange plotstop_row - cowplot::plot_grid(density_x, NULL,nrow 1, rel_widths c(3, 1), rel_heights c(0.5,0.5))bottom_row - cowplot::plot_grid(feature_scatter,density_y,nrow 1, rel_widths c(3, 1), rel_heights c(3,3))combined_plot - cowplot::plot_grid(top_row , bottom_row,nrow 2 , rel_widths c(3, 1), rel_heights c(0.25,3))ggsave(filename 20240221_Cebpd_Doublet.pdf, plot combined_plot, width 12, height 12, dpi 600)## Cilia-forming ##c(Pax8low/Prom1 Cilia-forming,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2)trans - subset(Epi_Named, idents c(Slc1a3med/Sox9 Cilia-forming))custom_labels - function(x) {ifelse(x %% 1 0, as.character(x), )
}feature_scatter - FeatureScatter( trans, Krt8,Col1a1,cols black) # Scatter plotNoLegend()labs(title NULL)theme(panel.grid.major element_line(color grey, size 0.5),panel.grid.minor element_blank())theme(axis.text.x element_text(color black, size 12),axis.title.x element_text(color black, size 14, face bold.italic),axis.ticks.x element_line(color black),axis.text.y element_text(color black, size 12),axis.title.y element_text(color black, size 14, face bold.italic),axis.ticks.y element_line(color black),plot.margin unit(c(0, 0, 0, 0), mm))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on x-axisscale_y_continuous(limits c(0,5) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on y-axisplot_data - as.data.frame(feature_scatter$data)# Create density plot for x-axis
density_x - ggplot(plot_data, aes(x Krt8 , fill black)) geom_density() theme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_blank(),axis.title.x element_blank(),axis.ticks.x element_blank(),axis.text.y element_text(color white, size 12),axis.title.y element_text(color white, size 14, face bold.italic),axis.ticks.y element_line(color white),plot.margin unit(c(0, 0, 0, 0), mm))scale_y_continuous(labels function(y) sprintf(%.0f, y)) NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on x-axis# Create density plot for y-axis
density_y - ggplot(plot_data, aes(x Col1a1 , fill black)) geom_density() coord_flip() # Flip axes for y-density plottheme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_text(color white, size 12),axis.title.x element_text(color white, size 14, face bold.italic),axis.ticks.x element_line(color white),axis.text.y element_blank(),axis.title.y element_blank(),axis.ticks.y element_blank(),plot.margin unit(c(0, 0, 0, 0), mm))NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,5))scale_y_continuous(limits c(0,1) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on y-axis# Arrange plotstop_row - cowplot::plot_grid(density_x, NULL,nrow 1, rel_widths c(3, 1), rel_heights c(0.5,0.5))bottom_row - cowplot::plot_grid(feature_scatter,density_y,nrow 1, rel_widths c(3, 1), rel_heights c(3,3))combined_plot - cowplot::plot_grid(top_row , bottom_row,nrow 2 , rel_widths c(3, 1), rel_heights c(0.25,3))ggsave(filename 20240221_SlcCilia_Doublet.pdf, plot combined_plot, width 12, height 12, dpi 600)## Pax8 Cilia-forming ##c(Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2)cancer - subset(Epi_Named, idents c(Pax8low/Prom1 Cilia-forming))custom_labels - function(x) {ifelse(x %% 1 0, as.character(x), )
}feature_scatter - FeatureScatter( cancer, Krt8,Col1a1,cols black) # Scatter plotNoLegend()labs(title NULL)theme(panel.grid.major element_line(color grey, size 0.5),panel.grid.minor element_blank())theme(axis.text.x element_text(color black, size 12),axis.title.x element_text(color black, size 14, face bold.italic),axis.ticks.x element_line(color black),axis.text.y element_text(color black, size 12),axis.title.y element_text(color black, size 14, face bold.italic),axis.ticks.y element_line(color black),plot.margin unit(c(0, 0, 0, 0), mm))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on x-axisscale_y_continuous(limits c(0,5) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on y-axisplot_data - as.data.frame(feature_scatter$data)# Create density plot for x-axis
density_x - ggplot(plot_data, aes(x Krt8 , fill black)) geom_density() theme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_blank(),axis.title.x element_blank(),axis.ticks.x element_blank(),axis.text.y element_text(color white, size 12),axis.title.y element_text(color white, size 14, face bold.italic),axis.ticks.y element_line(color white),plot.margin unit(c(0, 0, 0, 0), mm))scale_y_continuous(labels function(y) sprintf(%.0f, y)) NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on x-axis# Create density plot for y-axis
density_y - ggplot(plot_data, aes(x Col1a1 , fill black)) geom_density() coord_flip() # Flip axes for y-density plottheme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_text(color white, size 12),axis.title.x element_text(color white, size 14, face bold.italic),axis.ticks.x element_line(color white),axis.text.y element_blank(),axis.title.y element_blank(),axis.ticks.y element_blank(),plot.margin unit(c(0, 0, 0, 0), mm))NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,5))scale_y_continuous(limits c(0,1) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on y-axis# Arrange plotstop_row - cowplot::plot_grid(density_x, NULL,nrow 1, rel_widths c(3, 1), rel_heights c(0.5,0.5))bottom_row - cowplot::plot_grid(feature_scatter,density_y,nrow 1, rel_widths c(3, 1), rel_heights c(3,3))combined_plot - cowplot::plot_grid(top_row , bottom_row,nrow 2 , rel_widths c(3, 1), rel_heights c(0.25,3))ggsave(filename 20240221_Prom1_Doublet.pdf, plot combined_plot, width 12, height 12, dpi 600)## Spdef Secretory ##c(Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2)sec1 - subset(Epi_Named, idents c(Spdef Secretory))custom_labels - function(x) {ifelse(x %% 1 0, as.character(x), )
}feature_scatter - FeatureScatter( sec1, Krt8,Col1a1,cols black) # Scatter plotNoLegend()labs(title NULL)theme(panel.grid.major element_line(color grey, size 0.5),panel.grid.minor element_blank())theme(axis.text.x element_text(color black, size 12),axis.title.x element_text(color black, size 14, face bold.italic),axis.ticks.x element_line(color black),axis.text.y element_text(color black, size 12),axis.title.y element_text(color black, size 14, face bold.italic),axis.ticks.y element_line(color black),plot.margin unit(c(0, 0, 0, 0), mm))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on x-axisscale_y_continuous(limits c(0,5) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on y-axisplot_data - as.data.frame(feature_scatter$data)# Create density plot for x-axis
density_x - ggplot(plot_data, aes(x Krt8 , fill black)) geom_density() theme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_blank(),axis.title.x element_blank(),axis.ticks.x element_blank(),axis.text.y element_text(color white, size 12),axis.title.y element_text(color white, size 14, face bold.italic),axis.ticks.y element_line(color white),plot.margin unit(c(0, 0, 0, 0), mm))scale_y_continuous(labels function(y) sprintf(%.0f, y)) NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on x-axis# Create density plot for y-axis
density_y - ggplot(plot_data, aes(x Col1a1 , fill black)) geom_density() coord_flip() # Flip axes for y-density plottheme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_text(color white, size 12),axis.title.x element_text(color white, size 14, face bold.italic),axis.ticks.x element_line(color white),axis.text.y element_blank(),axis.title.y element_blank(),axis.ticks.y element_blank(),plot.margin unit(c(0, 0, 0, 0), mm))NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,5))scale_y_continuous(limits c(0,1) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on y-axis# Arrange plotstop_row - cowplot::plot_grid(density_x, NULL,nrow 1, rel_widths c(3, 1), rel_heights c(0.5,0.5))bottom_row - cowplot::plot_grid(feature_scatter,density_y,nrow 1, rel_widths c(3, 1), rel_heights c(3,3))combined_plot - cowplot::plot_grid(top_row , bottom_row,nrow 2 , rel_widths c(3, 1), rel_heights c(0.25,3))ggsave(filename 20240221_Spdef_Doublet.pdf, plot combined_plot, width 12, height 12, dpi 600)## Selenop Secretory ##c(Ciliated 1,Ciliated 2)sec2 - subset(Epi_Named, idents c(Selenop/Gstm2high Secretory))custom_labels - function(x) {ifelse(x %% 1 0, as.character(x), )
}feature_scatter - FeatureScatter( sec2, Krt8,Col1a1,cols black) # Scatter plotNoLegend()labs(title NULL)theme(panel.grid.major element_line(color grey, size 0.5),panel.grid.minor element_blank())theme(axis.text.x element_text(color black, size 12),axis.title.x element_text(color black, size 14, face bold.italic),axis.ticks.x element_line(color black),axis.text.y element_text(color black, size 12),axis.title.y element_text(color black, size 14, face bold.italic),axis.ticks.y element_line(color black),plot.margin unit(c(0, 0, 0, 0), mm))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on x-axisscale_y_continuous(limits c(0,5) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on y-axisplot_data - as.data.frame(feature_scatter$data)# Create density plot for x-axis
density_x - ggplot(plot_data, aes(x Krt8 , fill black)) geom_density() theme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_blank(),axis.title.x element_blank(),axis.ticks.x element_blank(),axis.text.y element_text(color white, size 12),axis.title.y element_text(color white, size 14, face bold.italic),axis.ticks.y element_line(color white),plot.margin unit(c(0, 0, 0, 0), mm))scale_y_continuous(labels function(y) sprintf(%.0f, y)) NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on x-axis# Create density plot for y-axis
density_y - ggplot(plot_data, aes(x Col1a1 , fill black)) geom_density() coord_flip() # Flip axes for y-density plottheme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_text(color white, size 12),axis.title.x element_text(color white, size 14, face bold.italic),axis.ticks.x element_line(color white),axis.text.y element_blank(),axis.title.y element_blank(),axis.ticks.y element_blank(),plot.margin unit(c(0, 0, 0, 0), mm))NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,5))scale_y_continuous(limits c(0,1) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on y-axis# Arrange plotstop_row - cowplot::plot_grid(density_x, NULL,nrow 1, rel_widths c(3, 1), rel_heights c(0.5,0.5))bottom_row - cowplot::plot_grid(feature_scatter,density_y,nrow 1, rel_widths c(3, 1), rel_heights c(3,3))combined_plot - cowplot::plot_grid(top_row , bottom_row,nrow 2 , rel_widths c(3, 1), rel_heights c(0.25,3))ggsave(filename 20240221_Selenop_Doublet.pdf, plot combined_plot, width 12, height 12, dpi 600)## Ciliated 1 ##c(Ciliated 2)cil1 - subset(Epi_Named, idents c(Ciliated 1))custom_labels - function(x) {ifelse(x %% 1 0, as.character(x), )
}feature_scatter - FeatureScatter( cil1, Krt8,Col1a1,cols black) # Scatter plotNoLegend()labs(title NULL)theme(panel.grid.major element_line(color grey, size 0.5),panel.grid.minor element_blank())theme(axis.text.x element_text(color black, size 12),axis.title.x element_text(color black, size 14, face bold.italic),axis.ticks.x element_line(color black),axis.text.y element_text(color black, size 12),axis.title.y element_text(color black, size 14, face bold.italic),axis.ticks.y element_line(color black),plot.margin unit(c(0, 0, 0, 0), mm))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on x-axisscale_y_continuous(limits c(0,5) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on y-axisplot_data - as.data.frame(feature_scatter$data)# Create density plot for x-axis
density_x - ggplot(plot_data, aes(x Krt8 , fill black)) geom_density() theme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_blank(),axis.title.x element_blank(),axis.ticks.x element_blank(),axis.text.y element_text(color white, size 12),axis.title.y element_text(color white, size 14, face bold.italic),axis.ticks.y element_line(color white),plot.margin unit(c(0, 0, 0, 0), mm))scale_y_continuous(labels function(y) sprintf(%.0f, y)) NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on x-axis# Create density plot for y-axis
density_y - ggplot(plot_data, aes(x Col1a1 , fill black)) geom_density() coord_flip() # Flip axes for y-density plottheme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_text(color white, size 12),axis.title.x element_text(color white, size 14, face bold.italic),axis.ticks.x element_line(color white),axis.text.y element_blank(),axis.title.y element_blank(),axis.ticks.y element_blank(),plot.margin unit(c(0, 0, 0, 0), mm))NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,5))scale_y_continuous(limits c(0,1) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on y-axis# Arrange plotstop_row - cowplot::plot_grid(density_x, NULL,nrow 1, rel_widths c(3, 1), rel_heights c(0.5,0.5))bottom_row - cowplot::plot_grid(feature_scatter,density_y,nrow 1, rel_widths c(3, 1), rel_heights c(3,3))combined_plot - cowplot::plot_grid(top_row , bottom_row,nrow 2 , rel_widths c(3, 1), rel_heights c(0.25,3))ggsave(filename 20240221_Ciliated_1_Doublet.pdf, plot combined_plot, width 12, height 12, dpi 600)## Ciliated 2 ##cil2 - subset(Epi_Named, idents c(Ciliated 2))custom_labels - function(x) {ifelse(x %% 1 0, as.character(x), )
}feature_scatter - FeatureScatter( cil2, Krt8,Col1a1,cols black) # Scatter plotNoLegend()labs(title NULL)theme(panel.grid.major element_line(color grey, size 0.5),panel.grid.minor element_blank())theme(axis.text.x element_text(color black, size 12),axis.title.x element_text(color black, size 14, face bold.italic),axis.ticks.x element_line(color black),axis.text.y element_text(color black, size 12),axis.title.y element_text(color black, size 14, face bold.italic),axis.ticks.y element_line(color black),plot.margin unit(c(0, 0, 0, 0), mm))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on x-axisscale_y_continuous(limits c(0,5) , breaks seq(0, 10, by 0.5),labels custom_labels) # Set breaks every 0.5 units on y-axisplot_data - as.data.frame(feature_scatter$data)# Create density plot for x-axis
density_x - ggplot(plot_data, aes(x Krt8 , fill black)) geom_density() theme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_blank(),axis.title.x element_blank(),axis.ticks.x element_blank(),axis.text.y element_text(color white, size 12),axis.title.y element_text(color white, size 14, face bold.italic),axis.ticks.y element_line(color white),plot.margin unit(c(0, 0, 0, 0), mm))scale_y_continuous(labels function(y) sprintf(%.0f, y)) NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,4) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on x-axis# Create density plot for y-axis
density_y - ggplot(plot_data, aes(x Col1a1 , fill black)) geom_density() coord_flip() # Flip axes for y-density plottheme(axis.line element_line(colorwhite),panel.background element_blank()) theme(axis.text.x element_text(color white, size 12),axis.title.x element_text(color white, size 14, face bold.italic),axis.ticks.x element_line(color white),axis.text.y element_blank(),axis.title.y element_blank(),axis.ticks.y element_blank(),plot.margin unit(c(0, 0, 0, 0), mm))NoLegend()scale_fill_manual(values c(grey))scale_x_continuous(limits c(0,5))scale_y_continuous(limits c(0,1) , breaks seq(0, 10, by 0.5)) # Set breaks every 0.5 units on y-axis# Arrange plotstop_row - cowplot::plot_grid(density_x, NULL,nrow 1, rel_widths c(3, 1), rel_heights c(0.5,0.5))bottom_row - cowplot::plot_grid(feature_scatter,density_y,nrow 1, rel_widths c(3, 1), rel_heights c(3,3))combined_plot - cowplot::plot_grid(top_row , bottom_row,nrow 2 , rel_widths c(3, 1), rel_heights c(0.25,3))ggsave(filename 20240221_Ciliated_2_Doublet.pdf, plot combined_plot, width 12, height 12, dpi 600)附图4 #### Figure Supp 4: Census of cell types of the mouse uterine tube ######## Packages Load ####library(dplyr)
library(patchwork)
library(Seurat)
library(harmony)
library(ggplot2)
library(cowplot)
library(SoupX)
library(DoubletFinder)
library(data.table)
library(parallel)
library(tidyverse)
library(SoupX)
library(ggrepel)library(ggplot2)
library(gplots)
library(RColorBrewer)
library(viridisLite)
library(Polychrome)
library(circlize)
library(NatParksPalettes)#### Proximal Datasets ####Proximal - readRDS( file ../dataset/Proximal_Filtered_Cells.rds , refhook NULL)Proximal_Named - RenameIdents(Proximal, 0 Fibroblast 1, 1 Stem-like Epithelial, 2 Fibroblast 2,3 Fibroblast 3, 4 Immune, 5 Secretory Epithelial, 6 Endothelial,7 Ciliated Epithelial,8 Mesothelial, 9 Smooth Muscle)Proximal_Namedactive.ident - factor(x Proximal_Namedactive.ident, levels c(Fibroblast 1,Fibroblast 2,Fibroblast 3,Smooth Muscle,Endothelial,Stem-like Epithelial,Secretory Epithelial,Ciliated Epithelial,Immune,Mesothelial))Proximal_Named - SetIdent(Proximal_Named, value Proximal_Namedactive.ident)Epi_Filter - readRDS(file ../dataset/Proximal_Epi_Cells.rds , refhook NULL)Epi_Named - RenameIdents(Epi_Filter, 0 Dbi/Spdefhigh Secretory, 1 Bmpr1b Progenitor, 2 Wfdc2 Secretory,3 Ciliated, 4 Sox17high Secretory, 5 Kcne3 Secretory)Epi_Namedactive.ident - factor(x Epi_Namedactive.ident, levels c(Ciliated,Dbi/Spdefhigh Secretory,Kcne3 Secretory,Sox17high Secretory,Wfdc2 Secretory,Bmpr1b Progenitor))#### Figure Supp 4a: Proximal All Cell Types ####Fibroblasts - c(#FF9D00 , #FFB653 , #FFCB9A) # Oranges
Muscle - c(#E55451) # Reds
Endothelial - c(#A0E6FF) # Blues
Epi -c(#6E3E6E,#CCCCFF,#DF73FF) # Purples
Immune - c( #5A5E6B ) # Grey
Meso - #1F51FF # Neon BLuecolors - c(Fibroblasts, Muscle, Endothelial, Epi, Immune, Meso)p1 - DimPlot(Proximal_Named,reductionumap,colscolors,pt.size 1.4,label.size 4,label.color black,repel TRUE,labelF) NoLegend() labs(xUMAP_1,yUMAP_2)ggsave(filename FIGs3a_all_proximal_umap.pdf, plot p1, width 15, height 12, dpi 600)#### Figure Supp 4b: Proximal Tile Mosaic ####library(treemap)Proximal_Named - RenameIdents(Proximal, 0 Fibroblast 1, 1 Stem-like Epithelial, 2 Fibroblast 2,3 Fibroblast 3, 4 Immune, 5 Secretory Epithelial, 6 Endothelial,7 Ciliated Epithelial,8 Mesothelial, 9 Smooth Muscle)Proximal_Namedactive.ident - factor(x Proximal_Namedactive.ident, levels c(Fibroblast 1,Fibroblast 2,Fibroblast 3,Smooth Muscle,Endothelial,Stem-like Epithelial,Secretory Epithelial,Ciliated Epithelial,Immune,Mesothelial))Proximal_Named - SetIdent(Proximal_Named, value Proximal_Namedactive.ident)Fibroblasts - c(#FF9D00 , #FFB653 , #FFCB9A) # Oranges
Muscle - c(#E55451) # Reds
Endothelial - c(#A0E6FF) # Blues
Epi -c(#6E3E6E,#CCCCFF,#DF73FF) # Purples
Immune - c( #5A5E6B ) # Grey
Meso - #1F51FF # Neon BLuecolors - c(Fibroblasts, Muscle, Endothelial, Epi, Immune, Meso)prox_cell_types - table(Idents(Proximal_Named), Proximal_Named$orig.ident)
prox_cell_type_df - as.data.frame(prox_cell_types)## Tile Mosaic ##prox_treemap - treemap(prox_cell_type_df, index Var1, vSize Freq, vColor colors, palette colors)ggsave(filename 20240612_all_prox_tile.pdf, plot prox_cell_type_df, width 8, height 12, dpi 600)#### Figure Supp 4c: Epi Markers All Distal Cells ####### Stacked Violin Plot Function ####https://divingintogeneticsandgenomics.rbind.io/post/stacked-violin-plot-for-visualizing-single-cell-data-in-seurat/## remove the x-axis text and tick
## plot.margin to adjust the white space between each plot.
## ... pass any arguments to VlnPlot in Seuratmodify_vlnplot - function(obj, feature, pt.size 0, plot.margin unit(c(-0.75, 0, -0.75, 0), cm),...) {p- VlnPlot(obj, features feature, pt.size pt.size, ... ) xlab() ylab(feature) ggtitle() theme(legend.position none, axis.text.x element_blank(), axis.ticks.x element_blank(), axis.title.y element_text(size rel(1), angle 0, face bold.italic), axis.text.y element_text(size rel(1)), plot.margin plot.margin ) return(p)
}## extract the max value of the y axis
extract_max- function(p){ymax- max(ggplot_build(p)$layout$panel_scales_y[[1]]$range$range)return(ceiling(ymax))
}## main function
StackedVlnPlot- function(obj, features,pt.size 0, plot.margin unit(c(-0.75, 0, -0.75, 0), cm),...) {plot_list- purrr::map(features, function(x) modify_vlnplot(obj obj,feature x, ...))# Add back x-axis title to bottom plot. patchwork is going to support this?plot_list[[length(plot_list)]]- plot_list[[length(plot_list)]] theme(axis.text.xelement_text(angle 60, hjust1, vjust0.95), axis.ticks.x element_line())# change the y-axis tick to only max value ymaxs- purrr::map_dbl(plot_list, extract_max)# plot_list- purrr::map2(plot_list, ymaxs, function(x,y) x plot_list- purrr::map2(plot_list, c(5,5,8,5), function(x,y) x scale_y_continuous(breaks c(y)) expand_limits(y y))p- patchwork::wrap_plots(plotlist plot_list, ncol 1)return(p)
}features- c(Epcam, Krt8 , Ovgp1 , Foxj1 )Fibroblasts - c(#FF9D00 , #FFB653 , #FFCB9A) # Oranges
Muscle - c(#E55451) # Reds
Endothelial - c(#A0E6FF) # Blues
Epi -c(#6E3E6E,#CCCCFF,#DF73FF) # Purples
Immune - c( #5A5E6B ) # Grey
Meso - #1F51FF # Neon BLuecolors - c(Fibroblasts, Muscle, Endothelial, Epi, Immune, Meso)stack_vln - StackedVlnPlot(obj Proximal_Named, features features, slot data,pt.size 0,cols colors) #9theme(plot.title element_text(size 32, face bold.italic))scale_x_discrete(limits c(Fibroblast 1,Fibroblast 2,Fibroblast 3,Smooth Muscle,Endothelial,Stem-like Epithelial,Secretory Epithelial,Ciliated Epithelial,Immune,Mesothelial))theme(axis.text.x element_text(size 16, angle 60))theme(axis.text.y element_text(size 14))theme(axis.title.y.left element_text(size 16))ggsave(filename 20240612_All_Prox_stacked_vln.pdf, plot stack_vln, width 18, height 12, dpi 600)#### Figure Supp 3d: Proximal Epi Cell Types ####epi_umap - DimPlot(object Epi_Named, # Seurat object reduction umap, # Axes for the plot (UMAP, PCA, etc.) #group.by Patient, # Labels to color the cells by (seurat_clusters, Age, Time.Point) repel TRUE, # Whether to repel the cluster labelslabel FALSE, # Whether to have cluster labels cols c( #35EFEF,#E95FE0,#B20224, #F28D86, #FB1111, #FEB0DB), pt.size 1.6, # Size of each dot is (0.1 is the smallest)label.size 2) # Font size for labels # You can add any ggplot2 1customizations herelabs(title Colored by Cluster) # Plot titleNoLegend()ggsave(filename FIGs3b_epi_proximal_umap.pdf, plot epi_umap, width 15, height 12, dpi 600)#### Figure Supp 4e: Epi Markers All Distal Cells ####P_Epi_Filter - readRDS(file ../Proximal/20220914_Proximal_Epi_Cells.rds , refhook NULL)P_Epi_Named - RenameIdents(P_Epi_Filter, 0 Dbi/Spdefhigh Secretory, 1 Bmpr1b Progenitor, 2 Wfdc2 Secretory,3 Ciliated, 4 Sox17high Secretory, 5 Kcne3 Secretory)P_Epi_Namedactive.ident - factor(x P_Epi_Namedactive.ident, levels c(Bmpr1b Progenitor,Wfdc2 Secretory,Sox17high Secretory,Kcne3 Secretory,Dbi/Spdefhigh Secretory,Ciliated))stack_vln - StackedVlnPlot(obj P_Epi_Named, features features, slot data,pt.size 0,cols c( #35EFEF,#F28D86, #FB1111,#FEB0DB,#B20224,#E95FE0))theme(plot.title element_text(size 32, face bold.italic))scale_x_discrete(limits c(Bmpr1b Progenitor,Wfdc2 Secretory,Sox17high Secretory,Kcne3 Secretory,Dbi/Spdefhigh Secretory,Ciliated))theme(axis.text.x element_text(size 16, angle 60))theme(axis.text.y element_text(size 14))theme(axis.title.y.left element_text(size 16))ggsave(filename 20240612_Prox_Epi_stacked_vln.pdf, plot stack_vln, width 18, height 12, dpi 600)#### Figure Supp 4f: Proximal Epi Features####named_features - c(Krt8,Epcam, Msln,Slc1a3,Sox9,Itga6, Bmpr1b,Ovgp1,Sox17,Pax8, Egr1,Wfdc2,Dbi,Gsto1,Fxyd4,Vim,Kcne3,Spdef,Lgals1,Upk1a, Thrsp,Selenop, Gstm2,Anpep, Klf6,Id2,Ifit1,Prom1, Ly6a, Kctd8, Adam8,Foxj1,Fam183b,Rgs22,Dnali1, Mt1 , Dynlrb2)prox_dp - DotPlot(object Epi_Named, # Seurat objectassay RNA, # Name of assay to use. Default is the active assayfeatures named_features, # List of features (select one from above or create a new one)# Colors to be used in the gradientcol.min 0, # Minimum scaled average expression threshold (everything smaller will be set to this)col.max 2.5, # Maximum scaled average expression threshold (everything larger will be set to this)dot.min 0, # The fraction of cells at which to draw the smallest dot (default is 0)dot.scale 6, # Scale the size of the pointsgroup.by NULL, # How the cells are going to be groupedsplit.by NULL, # Whether to split the data (if you fo this make sure you have a different color for each variable)scale TRUE, # Whether the data is scaledscale.by radius, # Scale the size of the points by size or radiusscale.min NA, # Set lower limit for scalingscale.max NA # Set upper limit for scaling
) labs(x NULL, # x-axis labely NULL)scale_color_viridis_c(optionF,begin.4,end0.9, direction -1)geom_point(aes(sizepct.exp), shape 21, colourblack, stroke0.6)theme_linedraw()guides(x guide_axis(angle 90)) theme(axis.text.x element_text(size 12 , face italic))theme(axis.text.y element_text(size 12))theme(legend.title element_text(size 12)) scale_y_discrete(limits c(Ciliated,Dbi/Spdefhigh Secretory,Kcne3 Secretory,Sox17high Secretory,Wfdc2 Secretory,Bmpr1b Progenitor))ggsave(filename FIGs3c_epi_proximal_dotplot.pdf, plot prox_dp, width 18, height 10, dpi 600)附图5
#### Figure Supp 5: Distal and proximal epithelial cell correlation ######## Packages Load ####
library(dplyr)
library(patchwork)
library(Seurat)
library(harmony)
library(ggplot2)
library(cowplot)
library(SoupX)
library(DoubletFinder)
library(data.table)
library(parallel)
library(tidyverse)
library(SoupX)
library(ggrepel)library(ggplot2)
library(gplots)
library(RColorBrewer)
library(viridisLite)
library(Polychrome)
library(circlize)
library(NatParksPalettes)#### Proximal Datasets ####Proximal - readRDS( file ../dataset/Proximal_Filtered_Cells.rds , refhook NULL)Proximal_Named - RenameIdents(Proximal, 0 Fibroblast 1, 1 Stem-like Epithelial, 2 Fibroblast 2,3 Fibroblast 3, 4 Immune, 5 Secretory Epithelial, 6 Endothelial,7 Ciliated Epithelial,8 Mesothelial, 9 Smooth Muscle)Proximal_Namedactive.ident - factor(x Proximal_Namedactive.ident, levels c(Fibroblast 1,Fibroblast 2,Fibroblast 3,Smooth Muscle,Endothelial,Stem-like Epithelial,Secretory Epithelial,Ciliated Epithelial,Immune,Mesothelial))Proximal_Named - SetIdent(Proximal_Named, value Proximal_Namedactive.ident)Epi_Filter - readRDS(file ../dataset/Proximal_Epi_Cells.rds , refhook NULL)Epi_Named - RenameIdents(Epi_Filter, 0 Dbi/Spdefhigh Secretory, 1 Bmpr1b Progenitor, 2 Wfdc2 Secretory,3 Ciliated, 4 Sox17high Secretory, 5 Kcne3 Secretory)Epi_Namedactive.ident - factor(x Epi_Namedactive.ident, levels c(Ciliated,Dbi/Spdefhigh Secretory,Kcne3 Secretory,Sox17high Secretory,Wfdc2 Secretory,Bmpr1b Progenitor))#### Proximal vs Distal Cluster Correlation ####Distal_Epi_Filter - readRDS(file ../dataset/Distal_Epi_Cells.rds , refhook NULL)Distal_Epi_Named - RenameIdents(Distal_Epi_Filter, 0 Spdef Secretory, 1 Slc1a3 Stem/Progenitor, 2 Cebpdhigh/Foxj1- Progenitor,3 Ciliated 1, 4 Ciliated 2, 5 Pax8low/Prom1 Cilia-forming, 6 Fibroblast-like,7 Slc1a3med/Sox9 Cilia-forming,8 Selenop/Gstm2high Secretory)Distal_Epi_Namedactive.ident - factor(x Distal_Epi_Namedactive.ident, levels c( c(Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming, Fibroblast-like,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2)))prox_avg_exp - AverageExpression(Epi_Named)$RNA
distal_avg_exp - AverageExpression(Distal_Epi_Named)$RNAcor.exp - as.data.frame(cor(x prox_avg_exp , y distal_avg_exp))cor.exp$x - rownames(cor.exp)cor.df - tidyr::gather(data cor.exp, y, correlation, c(Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming, Fibroblast-like,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2))distal_cells - c(Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming, Fibroblast-like,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2)prox_cells - c(Bmpr1b Progenitor,Ciliated,Dbi/Spdefhigh Secretory,Wfdc2 Secretory,Sox17high Secretory,Kcne3 Secretory)corr_matrix - ggplot(cor.df, aes(x, y, fill correlation)) geom_tile(color black,lwd 1,linetype 1) scale_fill_viridis_c(values c(0,1),optionrocket, begin.4,end0.99, direction -1,)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold), # Text size throughout the plotaxis.text.x element_text(color black, angle 60, hjust 1, size 12, face bold), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 12, face bold.italic))theme(plot.title element_blank())scale_y_discrete(limits c(Ciliated 2,Ciliated 1,Selenop/Gstm2high Secretory,Spdef Secretory,Fibroblast-like,Pax8low/Prom1 Cilia-forming, Slc1a3med/Sox9 Cilia-forming,Cebpdhigh/Foxj1- Progenitor,Slc1a3 Stem/Progenitor))scale_x_discrete(limits c(Bmpr1b Progenitor,Wfdc2 Secretory,Sox17high Secretory,Kcne3 Secretory,Dbi/Spdefhigh Secretory,Ciliated))geom_text(aes(x, y, label round(correlation, digits 2)), color black, size 4)ggsave(filename FIGs3d_epi_cluster_corr.pdf, plot corr_matrix, width 18, height 10, dpi 600)附图6 #### Figure Supp 6 and 10: Stem and Cancer Markers ######## Packages Load ####library(dplyr)
library(patchwork)
library(Seurat)
library(harmony)
library(ggplot2)
library(cowplot)
library(SoupX)
library(DoubletFinder)
library(data.table)
library(parallel)
library(tidyverse)
library(SoupX)
library(ggrepel)library(ggplot2)
library(gplots)
library(RColorBrewer)
library(viridisLite)
library(Polychrome)
library(circlize)
library(NatParksPalettes)library(monocle3)
library(ComplexHeatmap)
library(ggExtra)
library(gridExtra)
library(egg)library(scales)#### Distal Epithelial and Pseudotime Dataset ####Epi_Filter - readRDS(file ../dataset/Distal_Epi_Cells.rds , refhook NULL)Epi_Named - RenameIdents(Epi_Filter, 0 Spdef Secretory, 1 Slc1a3 Stem/Progenitor, 2 Cebpdhigh/Foxj1- Progenitor,3 Ciliated 1, 4 Ciliated 2, 5 Pax8low/Prom1 Cilia-forming, 6 Fibroblast-like,7 Slc1a3med/Sox9 Cilia-forming,8 Selenop/Gstm2high Secretory)Epi_Namedactive.ident - factor(x Epi_Namedactive.ident, levels c( c(Slc1a3 Stem/Progenitor,Cebpdhigh/Foxj1- Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming, Fibroblast-like,Spdef Secretory,Selenop/Gstm2high Secretory,Ciliated 1,Ciliated 2)))cds - readRDS(file ../dataset/Distal_Epi_PHATE_Monocle3.rds , refhook NULL)#### Figure Supp 4: Stem Dot Plot ####stem_features - c(Krt5,Krt17,Cd44,Prom1,Kit,Aldh1a1,Aldh1a2,Aldh1a3,Efnb1,Ephb1,Trp63,Sox2,Sox9,Klf4,Rnf43,Foxm1,Pax8,Nanog,Itga6,Psca,Tcf3,Tcf4,Nrp1,Slc1a3,Tnfrsf19,Smo,Lrig1,Ezh2,Egr1,Tacstd2,Dusp1,Slc38a2,Malat1,Btg2,Cdkn1c,Pdk4,Nedd9,Fos,Jun,Junb,Zfp36,Neat1,Gadd45g,Gadd45b)stem_dp - DotPlot(object Epi_Named, # Seurat objectassay RNA, # Name of assay to use. Default is the active assayfeatures stem_features, # List of features (select one from above or create a new one)# Colors to be used in the gradientcol.min 0, # Minimum scaled average expression threshold (everything smaller will be set to this)col.max 2.5, # Maximum scaled average expression threshold (everything larger will be set to this)dot.min 0, # The fraction of cells at which to draw the smallest dot (default is 0)dot.scale 6, # Scale the size of the pointsgroup.by NULL, # How the cells are going to be groupedsplit.by NULL, # Whether to split the data (if you fo this make sure you have a different color for each variable)scale TRUE, # Whether the data is scaledscale.by radius, # Scale the size of the points by size or radiusscale.min NA, # Set lower limit for scalingscale.max NA ) # Set upper limit for scalinglabs(x NULL, # x-axis labely NULL)scale_color_viridis_c(optionF,begin.4,end0.9, direction -1)geom_point(aes(sizepct.exp), shape 21, colourblack, stroke0.6)#theme_linedraw()guides(x guide_axis(angle 90))theme(axis.text.x element_text(size 8 , face italic))theme(axis.text.y element_text(size 9))theme(legend.title element_text(size 9))theme(legend.text element_text(size 8)) scale_y_discrete(limits c(Ciliated 2,Ciliated 1,Selenop/Gstm2high Secretory,Spdef Secretory,Fibroblast-like,Pax8low/Prom1 Cilia-forming, Slc1a3med/Sox9 Cilia-forming,Cebpdhigh/Foxj1- Progenitor,Slc1a3 Stem/Progenitor))ggsave(filename FIGs4_stem_dp.pdf, plot stem_dp, width 12, height 6, dpi 600)x - stem_dp$datawrite.csv( x , stem_dp_data.csv)#### Figure Supp 6: HGSC Driver Gene by Pseudotime ###### Calculate Pseudotime Values ##pseudo - pseudotime(cds)Distal_PHATEmeta.data$Pseudotime - pseudo # Add to Seurat Metadata## Subset Seurat Object ##color_cells - DimPlot(Distal_PHATE , reduction phate, cols c(#B20224, #1#35EFEF, #2#00A1C6, #3#A374B5, #4#9000C6, #5#EA68E1, #6lightgrey, #7#2188F7, #8#F28D86),pt.size 0.7,shuffle TRUE,seed 0,label FALSE)## Psuedotime and Lineage Assignment ##cellID - rownames(Distal_PHATEreductions$phatecell.embeddings)
phate_embeddings - Distal_PHATEreductions$phatecell.embeddings
pseudotime_vals - Distal_PHATEmeta.data$Pseudotimecombined_data - data.frame(cellID, phate_embeddings, pseudotime_vals)# Calculate the Average PHATE_1 Value for Pseudotime Points 0 #
avg_phate_1 - mean(phate_embeddings[pseudotime_vals 0, 1])# Pseudotime Values lower than avge PHATE_1 Embedding will be Negative to split lineages
combined_data$Split_Pseudo - ifelse(phate_embeddings[, 1] avg_phate_1, -pseudotime_vals, pseudotime_vals)# Define Lineage #
combined_data$lineage - ifelse(combined_data$PHATE_1 avg_phate_1, Secretory,ifelse(combined_data$PHATE_1 avg_phate_1, Ciliogenic, Progenitor))Distal_PHATE$Pseudotime_Adj - combined_data$Split_Pseudo
Distal_PHATE$Lineage - combined_data$lineage# Subset #Pseudotime_Lineage - subset(Distal_PHATE, idents c(Secretory 1,Secretory 2,Msln Progenitor,Slc1a3/Sox9 Cilia-forming,Pax8/Prom1 Cilia-forming,Progenitor,Ciliated 1,Ciliated 2))## Set Bins ##bins - cut_number(Pseudotime_Lineagemeta.data$Pseudotime_Adj , 40) # Evenly distribute bins Pseudotime_Lineagemeta.data$Bin - bins # Metadata for Bins## Set Idents to PSeudoime Bin ##time_ident - SetIdent(Pseudotime_Lineage, value Pseudotime_Lineagemeta.data$Bin)av.exp - AverageExpression(time_ident, return.seurat T)$RNA # Calculate Avg log normalized expression# Calculates Average Expression for Each Bin #
# if you set return.seuratT, NormalizeData is called which by default performs log-normalization #
# Reported as avg log normalized expression ### Pseudotime Scale Bar ##list - 1:40
colors c(rev(rainbow20),rainbow20)
df - data.frame(data list, color colors)pseudo_bar - ggplot(df, aes(x 1:40, y 1, fill color)) geom_bar(stat identity,position fill, color black, size 0, width 1) scale_fill_identity() theme_void() theme(axis.line element_blank(),axis.ticks element_blank(),axis.text element_blank(),axis.title element_blank())ggsave(filename pseudo_bar.pdf, plot pseudo_bar, width 0.98, height 0.19, dpi 600)## Plot HGSC driver gene list across pseudotime bin ##features - c(Trp53, Brca1, Brca2, Csmd3, Nf1, Fat3, Gabra6, Rb1, Apc, Lrp1b,Prim2, Cdkn2a, Crebbp, Wwox, Ankrd11, Map2k4, Fancm, Fancd2, Rad51c, Pten)# Create Bin List and expression of features #bin_list - unique(Pseudotime_Lineagemeta.data$Bin) plot_info - as.data.frame(av.exp[features,]) # Call Avg Expression for featuresz_score - transform(plot_info, SDapply(plot_info,1, mean, na.rm TRUE))
z_score - transform(z_score, MEANapply(plot_info,1, sd, na.rm TRUE))z_score1 - (plot_info-z_score$MEAN)/z_score$SDplot_info$y - rownames(plot_info) # y values as features
z_score1$y - rownames(plot_info)plot_info - gather(data plot_info, x, expression, bin_list) #set plot
z_score1 - gather(data z_score1, x, z_score, bin_list) #set plot# Create Cell Clusters DF #Labeled_Pseudotime_Lineage - RenameIdents(Pseudotime_Lineage, Secretory 1 Spdef Secretory, Progenitor Slc1a3 Stem/Progenitor, Msln Progenitor Cebpdhigh/Foxj1- Progenitor,Ciliated 1 Ciliated 1, Ciliated 2 Ciliated 2, Pax8/Prom1 Cilia-forming Pax8low/Prom1 Cilia-forming, Fibroblast-like Fibroblast-like, #removedSlc1a3/Sox9 Cilia-forming Slc1a3med/Sox9 Cilia-forming,Secretory 2 Selenop/Gstm2high Secretory)cluster_table - table(Labeled_Pseudotime_Lineageactive.ident, Labeled_Pseudotime_Lineagemeta.data$Bin)clusters - data.frame(cluster_table)clusters - clusters %% group_by(Var2) %%mutate(Perc Freq / sum(Freq))# Create Pseudotime DF #pseudotime_table - table(seq(1, length(bin_list), 1), unique(Labeled_Pseudotime_Lineagemeta.data$Bin),seq(1, length(bin_list), 1))pseudotime_bins - data.frame(pseudotime_table) # calculate max and min z-scores
max_z - max(z_score1$z_score, na.rm TRUE)
min_z - min(z_score1$z_score, na.rm TRUE)# set color for outliers
outlier_color - ifelse(z_score1$z_score max_z | z_score1$z_score min_z, ifelse(z_score1$z_score 0, #AD1F24, #51A6DC), #e2e2e2)## Plot Gene Expression ### Set different na.value options for positive and negative values
na_color_pos - #AD1F24 # color for positive NA values
na_color_neg - #51A6DC # color for negative NA valuescustom_bin_names - c(paste0(S, 20:1), paste0(C, 1:20))figure - ggplot(z_score1, aes(x, y, fill z_score)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradientn(colorsc(#1984c5, #e2e2e2, #c23728), name Average Expression \nZ-Score, limits c(-3,3), na.value ifelse(is.na(z_score1) z_score1 0, na_color_pos, ifelse(is.na(z_score1) z_score1 0, na_color_neg, grey50)),oob scales::squish)scale_x_discrete(limits sort(bin_list) , labels custom_bin_names)scale_y_discrete(limits rev(features))theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold), # Text size throughout the plotaxis.text.x element_text(color black, angle 0, hjust 0.5, size 10, face bold), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold.italic))theme(plot.title element_blank(),plot.marginunit(c(-0.5,1,1,1), cm))## Plot Cluster Percentage ##Spdef Secretory - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Spdef Secretory)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(1,1,1,1), cm))Selenop/Gstm2high Secretory - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Selenop/Gstm2high Secretory)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))Cebpdhigh/Foxj1- Progenitor - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Cebpdhigh/Foxj1- Progenitor)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))Slc1a3 Stem/Progenitor - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Slc1a3 Stem/Progenitor)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))Slc1a3med/Sox9 Cilia-forming - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Slc1a3med/Sox9 Cilia-forming)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))Pax8low/Prom1 Cilia-forming - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Pax8low/Prom1 Cilia-forming)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))Ciliated 1 - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Ciliated 1)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))Ciliated 2 - ggplot(clusters, aes(Var2, Var1, fill Perc)) geom_tile(color black,lwd 1,linetype 1) scale_fill_gradient2(lowwhite, high#000000, mid white, midpoint 0, name Percentage)scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Ciliated 2)theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))## Plot Pseudotime Color ##list - 1:40
colors c(rev(rainbow20),rainbow20)
df - data.frame(data list, color colors)binning - ggplot(df, aes(x 1:40, y 1, fill color)) geom_bar(stat identity,position fill, color black, size 1, width 1) scale_fill_identity() theme_void() theme(axis.line element_blank(),axis.ticks element_blank(),axis.text element_blank(),axis.title element_blank())scale_x_discrete(limits sort(bin_list) , labels seq(1, length(bin_list), 1))scale_y_discrete(limits Pseudotime Bin )theme(panel.background element_blank())labs(title Expression of Genes by Pseudotime Bin ,x element_blank(),y element_blank())theme(text element_text(size 12, face bold),# Text size throughout the plotaxis.ticks.x element_blank(),axis.ticks.y element_blank(),axis.text.x element_blank(), # Text color, angle, and horizontal adjustment on x-axis axis.text.y element_text(color black, hjust 1, vjust .75, size 14, face bold))theme(plot.title element_blank(),plot.marginunit(c(-1.25,1,1,1), cm))### Combine Plots ###psuedotime_lineage - ggarrange(Spdef Secretory,Selenop/Gstm2high Secretory,Cebpdhigh/Foxj1- Progenitor,Slc1a3 Stem/Progenitor,Slc1a3med/Sox9 Cilia-forming,Pax8low/Prom1 Cilia-forming,Ciliated 1,Ciliated 2,binning,figure , ncol1,heights c(2, 2, 2, 2, 2, 2, 2, 2, 2, (2*length(features)),widths c(3)),padding unit(0.01))ggsave(filename FIGs6_psuedotime_driver_gene.pdf, plot psuedotime_lineage, width 18, height 9, dpi 600)write.csv(z_score1 , cancer_pseudotime.csv)