Integrated Analysis Dotplot subscript out of bounds

[mpatino2799]

I am following the "Integrating stimulated vs. control PBMC datasets" tutorial for my own data comparing two conditions. All works well until I try to plot conserved markers with dotplot. I get the error:

Error in FUN(X[[i]], ...) : subscript out of bounds

Below is the code I tried to use for the dotplot:

Idents(RLU.combined) <- factor(Idents(RLU.combined), levels = c("1","8"))
markers.to.plot <- c("Gad1", "Sst")
DotPlot(
  RLU.combined,  features = rev(markers.to.plot), 
  cols = c("blue", "red"), dot.scale = 8, 
  split.by = "stim"
) + RotatedAxis()

Here is a link to the object I used: https://drive.google.com/open?id=13mlI9k1VFb8sLwoFNvwhi0U-H_YFih7Q

I downsampled the original object with RLU.combined.small <- subset(RLU.combined, downsample = 100)

[liu-xingliang]

By default, DotPlot use the data stored in data slot of current active assay.

I found your subset object has already with default assay set to RNA:

RLU.combined.small <- readRDS("RLU.combined.small.Rds")
# > dim(RLU.combined.small)
# Loading required package: Seurat
# [1] 25843    20
# > DefaultAssay(RLU.combined.small)
# [1] "RNA"

However, the gene name is in MGI format, not HGNC as you want. You have to update them.

# > RLU.combined.small@assays$RNA@data[1:3, 1:3]
# 3 x 3 sparse Matrix of class "dgCMatrix"
#               Calb2_tdTpositive_cell_46 PvalbD-Slc32a1_tdTpositive_cell_16
# 0610005C13Rik                 .                                  .        
# 0610007C21Rik                 0.8684538                          0.9306731
# 0610007L01Rik                 .                                  0.8672077
#               PvalbD-Slc32a1_tdTpositive_cell_31
# 0610005C13Rik                          .        
# 0610007C21Rik                          0.4667605
# 0610007L01Rik                          .        

You can do as below (take look at #978).

rownames(RLU.combined.small@assays$RNA@data) <- # new names

bless~

[andrewwbutler]

Hi,

The issue stems from one of your cells having a missing ident (NA), not a problem with the gene names. Removing (or reclassifying) this cell in the downsampled object resolves the problem. Please be sure every cell has been assigned a valid ident.

[mpatino2799]

Hi @andrewwbutler, thank you for taking a look. How do I remove or reclassify cells with missing identities?

Also for future reference, how can I check if cells have a missing ident (NA)?

I tried any(is.na(RLU.combined@meta.data[["orig.ident"]])) and it returns FALSE. I assume I am checking in the wrong place.

[liu-xingliang]

Sorry that I didn't check your data thoroughly. Thank u @andrewwbutler for pointing that out.

data <- readRDS("RLU.combined.small.Rds")
markers.to.plot <- c("Gad1", "Sst")

FetchData(data, vars=markers.to.plot)
#                                        Gad1        Sst
# PvalbD-Slc32a1_tdTpositive_cell_16 3.968287 0.00000000
# PvalbD-Slc32a1_tdTpositive_cell_31 3.955163 0.00260546
# Pvalb_tdTpositive_cell_25          3.118672 0.00000000
# PvalbF-Gad2_tdTpositive_cell_4     4.373908 0.00000000
# Pvalb_tdTpositive_cell_44          3.671748 0.00000000
# PvalbD-Slc32a1_tdTpositive_cell_39 3.634980 0.02060023
# Htr3a_tdTpositive_cell_85          4.157139 0.00000000
# Pvalb_tdTpositive_cell_64          3.412475 0.05999751
# PvalbF-Gad2_tdTpositive_cell_48    3.971190 0.00000000
# PvalbD-Slc32a1_tdTpositive_cell_12 3.273794 0.00000000
# Pvalb_tdTpositive_cell_5           3.438514 0.00000000
# PvalbF-Gad2_tdTpositive_cell_15    3.820229 0.00000000
# Sst_tdTpositive_cell_44            2.668655 3.42273627
# Sst_tdTpositive_cell_91            3.181188 2.93528798
# Sst_tdTpositive_cell_107           3.374119 3.14581586
# Sst_tdTpositive_cell_62            2.831219 3.85239623
# Htr3a_tdTpositive_cell_82          1.914746 3.60063334
# AGTGATCCATTGCTTT                   0.000000 0.00000000
# ATAGACCTCTCGACCT                   1.391159 0.00000000

# > Idents(data)
#          Calb2_tdTpositive_cell_46 PvalbD-Slc32a1_tdTpositive_cell_16 
#                               <NA>                                  1 
# PvalbD-Slc32a1_tdTpositive_cell_31          Pvalb_tdTpositive_cell_25 
#                                  1                                  1 
#     PvalbF-Gad2_tdTpositive_cell_4          Pvalb_tdTpositive_cell_44 
#                                  1                                  1 
# PvalbD-Slc32a1_tdTpositive_cell_39          Htr3a_tdTpositive_cell_85 
#                                  1                                  1 
#          Pvalb_tdTpositive_cell_64    PvalbF-Gad2_tdTpositive_cell_48 
#                                  1                                  1 
# PvalbD-Slc32a1_tdTpositive_cell_12           Pvalb_tdTpositive_cell_5 
#                                  1                                  1 
#    PvalbF-Gad2_tdTpositive_cell_15            Sst_tdTpositive_cell_44 
#                                  1                                  8 
#            Sst_tdTpositive_cell_91           Sst_tdTpositive_cell_107 
#                                  8                                  8 
#            Sst_tdTpositive_cell_62          Htr3a_tdTpositive_cell_82 
#                                  8                                  8 
#                   AGTGATCCATTGCTTT                   ATAGACCTCTCGACCT 
#                                  8                                  1 
# Levels: 1 8

data <- data[, !is.na(Idents(data))]

# > Idents(data)
# PvalbD-Slc32a1_tdTpositive_cell_16 PvalbD-Slc32a1_tdTpositive_cell_31 
#                                  1                                  1 
#          Pvalb_tdTpositive_cell_25     PvalbF-Gad2_tdTpositive_cell_4 
#                                  1                                  1 
#          Pvalb_tdTpositive_cell_44 PvalbD-Slc32a1_tdTpositive_cell_39 
#                                  1                                  1 
#          Htr3a_tdTpositive_cell_85          Pvalb_tdTpositive_cell_64 
#                                  1                                  1 
#    PvalbF-Gad2_tdTpositive_cell_48 PvalbD-Slc32a1_tdTpositive_cell_12 
#                                  1                                  1 
#           Pvalb_tdTpositive_cell_5    PvalbF-Gad2_tdTpositive_cell_15 
#                                  1                                  1 
#            Sst_tdTpositive_cell_44            Sst_tdTpositive_cell_91 
#                                  8                                  8 
#           Sst_tdTpositive_cell_107            Sst_tdTpositive_cell_62 
#                                  8                                  8 
#          Htr3a_tdTpositive_cell_82                   AGTGATCCATTGCTTT 
#                                  8                                  8 
#                   ATAGACCTCTCGACCT 
#                                  1 
# Levels: 1 8

DotPlot(
  data,  features = rev(markers.to.plot), 
  cols = c("blue", "red"), dot.scale = 8, 
  split.by = "stim"
) + RotatedAxis()