haojiawu / plot1cell

A R package for advanced single cell data viz

plot1cell: a package for advanced single cell data visualization

This R package allows users to visualize the single cell data on the R object or output files generated by Seurat. It is currently under active development.

Installation

plot1cell R package can be easily installed from Github using devtools. Please make sure you have installed Seurat 4.0, circlize and ComplexHeatmap packages.

devtools::install_github("TheHumphreysLab/plot1cell")
## or the development version, devtools::install_github("HaojiaWu/plot1cell")

## You might need to install the dependencies below if they are not available in your R library.
bioc.packages <- c("biomaRt","GenomeInfoDb","EnsDb.Hsapiens.v86","GEOquery","simplifyEnrichment","ComplexHeatmap")
BiocManager::install(bioc.packages)
dev.packages <- c("chris-mcginnis-ucsf/DoubletFinder","Novartis/hdf5r","mojaveazure/loomR")
devtools::install_github(dev.packages)
## If you can't get the hdf5r package installed, please see the fix here:
## https://github.com/hhoeflin/hdf5r/issues/94

Usage

We provide some example codes to help generate figures from user's provided Seurat object. The Seurat object input to plot1cell should be a final object with complete clustering and cell type annotation. If a seurat object is not available, we suggest to use the demo data from Satija's lab (https://satijalab.org/seurat/articles/integration_introduction.html). To demonstrate the plotting functions in plot1cell, we re-created a Seurat object from our recent paper Kirita et al, PNAS 2020 by integrating the count matrices we uploaded to GEO (GSE139107).

library(plot1cell)
iri.integrated <- Install.example() 

# Please note that this Seurat object is just for demo purpose and 
# is not exactly the same as the one we published on PNAS.
# It takes about 2 hours to run in a linux server with 500GB RAM and 32 CPU cores.
# You can skip this step and use your own Seurat object instead

1. Circlize plot to visualize cell clustering and meta data

This circlize plot was inspired by the data visualization in a published paper (Figure1, https://www.nature.com/articles/s41586-021-03775-x) from Linnarsson's lab.

###Check and see the meta data info on your Seurat object
colnames([e-mail])  

###Prepare data for ploting
circ_data <- prepare_circlize_data(iri.integrated, scale = 0.8 )
set.seed(1234)
cluster_colors<-rand_color(length(levels(iri.integrated)))
group_colors<-rand_color(length(names(table(iri.integrated$Group))))
rep_colors<-rand_color(length(names(table(iri.integrated$orig.ident))))

###plot and save figures
png(filename =  'circlize_plot.png', width = 6, height = 6,units = 'in', res = 300)
plot_circlize(circ_data,do.label = T, pt.size = 0.01, col.use = cluster_colors ,bg.color = 'white', kde2d.n = 200, repel = T, label.cex = 0.6)
add_track(circ_data, group = "Group", colors = group_colors, track_num = 2) ## can change it to one of the columns in the meta data of your seurat object
add_track(circ_data, group = "orig.ident",colors = rep_colors, track_num = 3) ## can change it to one of the columns in the meta data of your seurat object
dev.off()

2. Dotplot to show gene expression across groups

Here is an example to use plot1cell to show one gene expression across different cell types in different groups.

png(filename =  'dotplot_single.png', width = 4, height = 6,units = 'in', res = 100)
complex_dotplot_single(seu_obj = iri.integrated, feature = "Havcr1",groups = "Group")
dev.off()

If the group factor can be classified by another factor, complex_dotplot_single allows splitting the group factor by another group factor too. Here is an example for demo.

[e-mail]$Phase<-plyr::mapvalues([e-mail]$Group, from = levels([e-mail]$Group), to = c("Healthy",rep("Injury",3), rep("Recovery",2)))
[e-mail]$Phase<-as.character([e-mail]$Phase)
png(filename =  'dotplot_single_split.png', width = 4, height = 6,units = 'in', res = 100)
complex_dotplot_single(iri.integrated, feature = "Havcr1",groups = "Group",splitby = "Phase")
dev.off()

To visualize the same gene on multiple group factors, simply add more group factor IDs to the groups argument.

png(filename =  'dotplot_more_groups.png', width = 8, height = 6,units = 'in', res = 100)
complex_dotplot_single(seu_obj = iri.integrated, feature = "Havcr1",groups= c("Group","Replicates"))
dev.off()

Each group factor can be further splitted by its own factor if the splitby argument is provided. Note that in this case, the order of the group factors needs to match the order of splitby factors.

[e-mail]$ReplicateID<-plyr::mapvalues([e-mail]$Replicates, from = names(table(([e-mail]$Replicates))), to = c(rep("Rep1",3),rep("Rep2",3), rep("Rep3",1)))
[e-mail]$ReplicateID<-as.character([e-mail]$ReplicateID)

png(filename =  'dotplot_more_groups_split.png', width = 9, height = 6,units = 'in', res = 200)
complex_dotplot_single(seu_obj = iri.integrated, feature = "Havcr1",groups= c("Group","Replicates"), splitby = c("Phase","ReplicateID"))
dev.off()
### In this example, "Phase" is a splitby factor for "Group" and "ReplicateID" is a splitby factor for "Replicates".

Note that the Replicates group here is just for showcase purpose. This is not a meaningful group ID in our snRNA-seq dataset.

To visualize multiple genes in dotplot format, complex_dotplot_multiple should be used.

png(filename =  'dotplot_multiple.png', width = 10, height = 4,units = 'in', res = 300)
complex_dotplot_multiple(seu_obj = iri.integrated, features = c("Slc34a1","Slc7a13","Havcr1","Krt20","Vcam1"),group = "Group", celltypes = c("PTS1" ,   "PTS2"  ,  "PTS3"  ,  "NewPT1" , "NewPT2"))
dev.off()

3. Violin plot to show gene expression across groups

One gene/one group factor violin plot:

png(filename =  'vlnplot_single.png', width = 4, height = 6,units = 'in', res = 100)
complex_vlnplot_single(iri.integrated, feature = "Havcr1", groups = "Group",celltypes   = c("PTS1" ,   "PTS2"  ,  "PTS3"  ,  "NewPT1" , "NewPT2"))
dev.off()

Similar to complex_dotplot_single, the complex_vlnplot_single function also allows splitting the group factor by another factor with the argument splitby.

png(filename =  'vlnplot_single_split.png', width = 4, height = 6,units = 'in', res = 100)
complex_vlnplot_single(iri.integrated, feature = "Havcr1", groups = "Group",celltypes   = c("PTS1" ,   "PTS2"  ,  "PTS3"  ,  "NewPT1" , "NewPT2"), splitby = "Phase")
dev.off()

One gene/multiple group factors violin plot:

png(filename =  'vlnplot_multiple.png', width = 6, height = 6,units = 'in', res = 100)
complex_vlnplot_single(iri.integrated, feature = "Havcr1", groups = c("Group","Replicates"),celltypes   = c("PTS1" ,   "PTS2"  ,  "PTS3"  ,  "NewPT1" , "NewPT2"), font.size = 10)
dev.off()

Similar to the functionality in complex_dotplot, each group factor can also be splitted by another factor in violin plot. For example:

png(filename =  'vlnplot_multiple_split.png', width = 7, height = 5,units = 'in', res = 200)
complex_vlnplot_single(iri.integrated, feature = "Havcr1", groups = c("Group","Replicates"),
                        celltypes   = c("PTS1" ,   "PTS2"  ,  "PTS3"  ,  "NewPT1" , "NewPT2"), 
                        font.size = 10, splitby = c("Phase","ReplicateID"), pt.size=0.05)
dev.off()

Multiple genes/one group factor violin plot:

png(filename =  'vlnplot_multiple_genes.png', width = 6, height = 6,units = 'in', res = 300)
complex_vlnplot_multiple(iri.integrated, features = c("Havcr1",  "Slc34a1", "Vcam1",   "Krt20"  , "Slc7a13", "Slc5a12"), celltypes = c("PTS1" ,   "PTS2"  ,  "PTS3"  ,  "NewPT1" , "NewPT2"), group = "Group", add.dot=T, pt.size=0.01, alpha=0.01, font.size = 10)
dev.off()

Multiple genes/multiple group factors.

The violin plot will look too messy in this scenario so it is not included in plot1cell.

4. Umap geneplot across groups

png(filename =  'data/geneplot_umap.png', width = 8, height = 6,units = 'in', res = 100)
complex_featureplot(iri.integrated, features = c("Havcr1",  "Slc34a1", "Vcam1",   "Krt20"  , "Slc7a13"), group = "Group", select = c("Control","12hours","6weeks"), order = F)
dev.off()

5. ComplexHeatmap to show unique genes across groups

plot1cell can directly identify the condition specific genes in a selected cell type and plot those genes using ComplexHeatmap. An example is shown below:

iri.integrated$Group2<-plyr::mapvalues(iri.integrated$Group, from = c("Control", "4hours",  "12hours", "2days",   "14days" , "6weeks" ),
to = c("Ctrl","Hr4","Hr12","Day2", "Day14","Wk6"))
iri.integrated$Group2<-factor(iri.integrated$Group2, levels = c("Ctrl","Hr4","Hr12","Day2", "Day14","Wk6"))
png(filename =  'heatmap_group.png', width = 4, height = 8,units = 'in', res = 100)
complex_heatmap_unique(seu_obj = iri.integrated, celltype = "NewPT2", group = "Group2",gene_highlight = c("Slc22a28","Vcam1","Krt20","Havcr1"))
dev.off()

6. Upset plot to show the unique and shared DEGs across groups.

png(filename =  'upset_plot.png', width = 8, height = 4,units = 'in', res = 300)
complex_upset_plot(iri.integrated, celltype = "NewPT2", group = "Group", min_size = 10, logfc=0.5)
dev.off()

7. Cell proportion change across groups

png(filename =  'cell_fraction.png', width = 8, height = 4,units = 'in', res = 300)
plot_cell_fraction(iri.integrated,  celltypes = c("PTS1" ,   "PTS2"  ,  "PTS3"  ,  "NewPT1" , "NewPT2"), groupby = "Group", show_replicate = T, rep_colname = "orig.ident")
dev.off()

8. Other ploting functions

There are other functions for plotting/data processing in plot1cell.

help(package = plot1cell)

Many more functions will be added in the future package development. For questions, please raise an issue in this github page or contact TheHumphreysLab.

9. Attributions

This package uses many methods from Seurat (https://github.com/satijalab/seurat) to process the data for ploting. The circlize and heatmap plots were generated by the circlize (https://github.com/jokergoo/circlize) and ComplexHeatmap (https://github.com/jokergoo/ComplexHeatmap) packages. The Upset plot was generated by the ComplexUpset package (https://github.com/krassowski/complex-upset). Most of other graphs were generated using ggplot2 (https://github.com/tidyverse/ggplot2). The package benefits from the following dependencies.

    Seurat,
    plotly,
    circlize,
    dplyr,
    ggplot2,
    ggh4x,
    MASS,
    scales,
    progress,
    RColorBrewer,
    grid,
    grDevices,
    biomaRt,
    reshape2,
    ggbeeswarm,
    purrr,
    ComplexUpset,
    matrixStats,
    DoubletFinder,
    methods,
    data.table,
    Matrix,
    hdf5r,
    loomR,
    GenomeInfoDb,
    EnsDb.Hsapiens.v86,
    cowplot,
    rlang,
    GEOquery,
    simplifyEnrichment,
    wordcloud,
    ComplexHeatmap

10. Citation

Please consider citing our paper if you find plot1cell useful.
https://www.cell.com/cell-metabolism/fulltext/S1550-4131(22)00192-9
Cell Metab. 2022 Jul 5;34(7):1064-1078.e6.
Wu H, Gonzalez Villalobos R, Yao X, Reilly D, Chen T, Rankin M, Myshkin E, Breyer MD, Humphreys BD.
Mapping the single-cell transcriptomic response of murine diabetic kidney disease to therapies.

plot1cell's Issues

"colors" setting in "add_track"

Thank you for the development of Plot1Cell. This issue is related to issue #6.

Here, I set the collor for plot_circlize and add_track. The plot_circlize works perfectly. But the add_track only use the ordered colors for each clusters, rather than use the correct color for each subsets clusters.

The code is as the following :
`
Order_Celltype_Main <- c(
"EC","VSMC",
"FC","CPC","PreHTC","RegC",
"M1 macrophage"
)
color_Celltype_Main <- c(
"#F8766D","#C49A00",
"#53B400","#00C094","#00B6EB","#A58AFF",
"#FB61D7"
)
names(color_Celltype_Main) <- Order_Celltype_Main
Seurat.object$Celltype_Main <- factor(
x = Seurat.object$Celltype_Main,
levels = Order_Celltype_Main
)

Order_Celltype_Subset <- c(
"EC (SELE)",
"VSMC (IL6)","VSMC (RGS5)",

"FC (AQP3)","FC (CFH)",
"FC (ECRG4)","FC (FOS)",
"FC (HSPB1)","FC (POSTN)",
"FC (PRG4)","FC (THBS4)",
"FC (LOC119877775)","FC (MAGI2)",
"FC (MMP3)","FC (MT1E)",
"FC (RBP4)","FC (S100A12)",

"CPC (STMN1)",
"PreHTC (MMP3)",
"RegC (PTGS2)",

"M1 macrophage (CCL23)","M1 macrophage"
)
color_Celltype_Subset <- c(
"#F8766D",
"#e0b311", "#b4900e",

"#9DB400", "#90B405", "#84B40B", "#78B411",
"#6CB417", "#60B41D", "#54B423", "#48B428", "#3CB42E",
"#30B434", "#24B43A", "#18B440", "#0CB446", "#00B44C",

"#00C094",
"#00B6EB",
"#A58AFF",

"#dd40b8", "#df74c6"
)
names(color_Celltype_Subset) <- Order_Celltype_Subset
Seurat.object$Celltype_Subset <- factor(
x = Seurat.object$Celltype_Subset,
levels = Order_Celltype_Subset
)

Plot1Cell.object <- prepare_circlize_data(Seurat.object, scale = 0.8 )

Plot1Cell.object$Celltype_Main <- factor(
x = Plot1Cell.object$Celltype_Main,
levels = Order_Celltype_Main
)
Plot1Cell.object$Celltype_Subset <- factor(
x = Plot1Cell.object$Celltype_Subset,
levels = Order_Celltype_Subset
)
Plot1Cell.object$Sample <- factor(
x = Plot1Cell.object$Sample,
levels = Order_Sample
)

plot_circlize(
Plot1Cell.object,
do.label = T,
contour.levels = c(0.1, 0.5, 0.9),
pt.size = 0.25,
col.use = color_Celltype_Main ,
bg.color = 'white',
kde2d.n = 1000,
repel = T,
label.cex = 0.7)

add_track(
Plot1Cell.object,
group = "Celltype_Subset",
colors = color_Celltype_Subset,
track_num = 2)
`

Error in the expression

Hi,

Thank you for the wonderful package. I was plotting the expression of genes using the fucntion:
complex_dotplot_single(seu_obj =sobj, feature = "gene1", groups = "condition")

I get no expression of gene1 in the last column. However, when I looked at the feature plot, I get a very high expression of the gene1 in the clusters. What do you think could be happening here?

Thank you

Error in anno_zoom using complex_heatmap_unique()

I encountered an error when I used "complex_heatmap_unique" function as follows:

complex_heatmap_unique(sce4, celltype = "10", group = "Group", gene_highlight = c("CDK1"))
Error in anno_zoom(...) :
unused arguments (add_new_line = TRUE, value_range = c(1, 2))

Here are the versions of the main packages I am using. Could you please take a look and see what might be the reason?
Thanks you so much

plot1cell_0.0.0.9000 ComplexHeatmap_2.10.0 Seurat_5.0.1 SeuratObject_5.0.1

Error while using the complex_dotplot_single() function

Hi @HaojiaWu,

Thank you so much for developing this package with the new plotting capabilities. I am trying to follow the tutorial and running into an error. I have a very large seurat object so I subsetted the object before trying to plot. My code is pasted below with the error message:

#Saving the plot since I keep getting an error in RStudio to plot in the plotting window
ggsave(complex_dotplot_single(cd8, feature = "HAVCR2",groups = [email protected]$Treatment,splitby = "Timepoint"),"test.png",device = "png")

##Error is pasted below
Warning: Found the following features in more than one assay, excluding the default. We will not include these in the final data frame: HAVCR2
Warning: Could not find RELA in the default search locations, found in RNA assay instead
Error in names(x) <- value :
'names' attribute [36604] must be the same length as the vector [3]
In addition: Warning message:
In FetchData.Seurat(new_seu, slot = data.type, vars = c(features, :
The following requested variables were not found: HAVCR2, NIVO, IPI

Please let me know if I am missing something that is triggering this error. If I use the same command without the "splitby" argument, it works just fine, and I am able to generate the plot.

Thanks a lot.

add_track width etc

Hi, I wonder if any parameter has been implemented to modify the width, and font size in the add_track function? Thanks.

Error: Cannot find 'umap' in this Seurat object

Hi!
This is my code
complex_dotplot_single(seu_obj = seurat_object,
feature = "Havcr1",
groups = "group2")
Error: Cannot find 'umap' in this Seurat object
Dimensional reduction for visualization was performed using tSNE in my seurat_object. Could you help solve this problem.
Thank you!

UMAP contours get cut off

Hi,
Thanks for developing this great visualization tool. When plotting the circlize plot, I am noticing that the UMAP contours are getting cut off in a square shape within the circle. Here are some examples from my plot.

I have tried adjusting the scale value in prepare_circlize_data to make it smaller, but the clipped areas remain. I have also noticed that your example on the front page also has the clipped edges. Any suggestions on how to fix this?

Thanks!

object 'UMAP1' not found when running complex_featureplot

When I ran the complex_featureplot function, I encountered the following error:

complex_featureplot(object, features = genes, group = "Group", order = T)

"Error in geom_point():
! Problem while computing aesthetics.
i Error occurred in the 1st layer.
Caused by error in FUN():
! object 'UMAP1' not found
Run rlang::last_trace() to see where the error occurred."

However, my object@reductions contains UMAP information.

I used the following code to change 'umap' to 'UMAP', yet the error persists.
object@reductions$UMAP <- object@reductions$umap

How can I resolve this issue? Thank you.

UMAP overlap with circ

Any parameter can adjust circle size?