Cryptic error on VlnPlot invalid 'times' argument

[rbutleriii]

I am trying to print a multipage pdf of stacked violin plots, for several gene sets. It normally works when all genes are found in variable count. However, in some cases it errors out and I would like to keep going. Toy example below

tasic_cort = c("Mybpc1", "Parm1", "Sncg", "Chat", "Gpc3", "Car4", "Cxcl14", 
               "Igtp", "Smad3", "Sncg", "Cbln4", "Th", "Myh8", "Cdk6", 
               "Tacstd2", "Chodl", "Gpx3", "Tpbg", "Cpne5", "Tacr3", "Rspo2", 
               "Wt1", "Obox3", "Ctxn3", "Scnn1a", "Hsd11b1", "Arf5", "Ptgs2", 
               "Ngb", "Tcerg1l", "Pde1c", "Batf3", "Ucma", "Car12", "Syt17", 
               "Rgs12", "Mgp", "Sla", "Chrna6", "Tph2", "Cdh13", 
               "Aqp4", "Pdgfra", "9630013A20Rik", "Opalin", "Ctss", "Xdh", 
               "Myl9")
DefaultAssay(sc) = "SCT"
CairoPDF(file=paste(today, nameset, type, "markers.pdf", sep="."), onefile=T, 
         width=11, height=8.5)
VlnPlot(sc, features=tasic_cort, stack=T) + 
        theme(legend.position="none") + ggtitle("Tasic etal 2016 All")
dev.off()

I think issue is a gene that is found with no counts (in the RNA_assay). As the errors look like this:

Warning: Could not find Tacstd2 in the default search locations, found in RNA assay instead
Warning: Could not find Wt1 in the default search locations, found in RNA assay instead
Warning: Could not find Obox3 in the default search locations, found in RNA assay instead
Warning: Could not find Chrna6 in the default search locations, found in RNA assay instead
Warning: All cells have the same value of rna_Wt1
Error in rep(x = "", times = length(x) - 2) : invalid 'times' argument
> traceback()
19: f(...)
18: self$labels(breaks)
17: f(..., self = self)
16: self$scale$get_labels(breaks)
15: f(..., self = self)
14: scale$get_labels(breaks)
13: guide_train.axis(guide, panel_params[[aesthetic]])
12: guide_train(guide, panel_params[[aesthetic]])
11: FUN(X[[i]], ...)
10: lapply(aesthetics, function(aesthetic) {
        axis <- substr(aesthetic, 1, 1)
        guide <- panel_params$guides[[aesthetic]]
        guide <- guide_train(guide, panel_params[[aesthetic]])
        guide <- guide_transform(guide, self, panel_params)
        guide <- guide_geom(guide, layers, default_mapping)
        guide
    })
9: f(..., self = self)
8: FUN(X[[i]], ...)
7: lapply(self$panel_params, self$coord$train_panel_guides, layers,
       default_mapping, self$coord_params)
6: f(..., self = self)
5: layout$setup_panel_guides(plot$guides, plot$layers, plot$mapping)
4: ggplot_gtable.ggplot_built(data)
3: ggplot_gtable(data)
2: print.ggplot(x)
1: (function (x, ...)
   UseMethod("print"))(x)

Maybe it cannot create axis breaks when all values are 0? Is there a way to tell VlnPlot to ignore missing gene names?

[longmanz]

Hi @rbutleriii, You may try running AverageExpression() first to filter the genes. For example:

ave_exp = AverageExpression(object = sc, assays = "RNA", features = tasic_cort)

You can remove the genes with 0 expression level (i.e., 0 count), then proceed to plotting. Best, Longda

[rbutleriii]

The simplest workaround I have found was just to look at the intersection of the gene list and the rownames of the active assay. This in effect disables looking at assay="RNA" and it simply excludes the genes missing in SCT:

VlnPlot(sc, features=intersect(tasic_cort, rownames(sc)), stack=T) + 
        theme(legend.position="none") + ggtitle("Tasic etal 2016 All")

But that obviously won't work if someone is plotting the RNA assay. Adding extra steps to the workflow to filter prior to each violin plot is ok...but isn't very practical if trying to use VlnPlot programmatically. In this case, I have an overall SeuratObject (SO) that I then subset by broad class and re-cluster. I can be virtually guaranteed the subsets will have zero count genes, but dropping them from the child SOs does become a minor annoyance when transferring annotations from the subset back to the main dataset.

Ideally, since the VlnPlot function has that multi-assay versatility built into it, it would be good if it didn't return an error and stop the Rscript, but just warned and dropped the gene from the plot.

[Sophia409]

I met the same question. Could the authors modify the code to make the VlnPlot function warn and drop the gene from the plot, not output an error?

[rbutleriii]

To touch back on this since I was just plotting the RNA assay. The easiest way I have found to exclude the zero count is just rowSums(sc) > 0. So to run through a stack of gene lists for comparison:

sc_names = rownames(sc)[rowSums(sc) > 0] 
CairoPDF(file=paste(today, nameset, typeout, "markers.pdf", sep="."), onefile=T, 
         width=11, height=8.5)
VlnPlot(sc, features=intersect(yang_cort, sc_names), stack=T, flip=T, raster=T) + 
        theme(legend.position="none") + ggtitle("Yang etal 2021 Cortex")
VlnPlot(sc, features=intersect(yang_chorplex, sc_names), stack=T, flip=T, raster=T) + 
        theme(legend.position="none") + ggtitle("Yang etal 2021 Choroid Plexus")
VlnPlot(sc, features=intersect(ximerakis, sc_names), stack=T, flip=T, raster=T) + 
        theme(legend.position="none") + ggtitle("Ximerakis etal 2019")
VlnPlot(sc, features=intersect(lake_all, sc_names), stack=T, raster=T) + 
        theme(legend.position="none") + ggtitle("Lake etal 2018 All")
VlnPlot(sc, features=intersect(lake_ex, sc_names), stack=T, raster=T) + 
        theme(legend.position="none") + ggtitle("Lake etal 2018 Excitatory Neurons")
VlnPlot(sc, features=intersect(lake_int, sc_names), stack=T, raster=T) + 
        theme(legend.position="none") + ggtitle("Lake etal 2018 Interneurons")
VlnPlot(sc, features=intersect(tasic_cort, sc_names), stack=T, raster=T) + 
        theme(legend.position="none") + ggtitle("Tasic etal 2016 All")
VlnPlot(sc, features=intersect(allen_cort, sc_names), stack=T, flip=T, raster=T) + 
        theme(legend.position="none") + ggtitle("Allen Brain Atlas 2019 Cortex")
dev.off() # filesize is still huge, look at lowering resolution

[SebastianMHJohn]

For me the scCustomize version worked fine: https://samuel-marsh.github.io/scCustomize/reference/Stacked_VlnPlot.html

Best wishes