Cannot RUN PCA using the codes of Cell_Plex_10X_Genomics with Seurat R

Question

[shweta.johari13]

Hello All, I am running a CMO data using the resources of 10X Genomics and combining it with Seurat R for further analysis as mentioned by the 10X Genomics. I am receiving an error and cannot perform RUNPCA command. Can anyone help me please? Please find enclosed the 10X Genomics code: https://www.10xgenomics.com/resources/analysis-guides/tag-assignment-of-10x-Genomics-cellplex-data-using-seurats-htodemux-function They have attached the vignette of Seurat: https://satijalab.org/seurat/articles/hashing_vignette.html

My entire Code:

setwd("E:/Cell_Ranger_10x_Seurat")
data_dir <- "E:/Cell_Ranger_10x_Seurat"
data <- Read10X(data.dir = data_dir)
10X data contains more than one type and is being returned as a list containing matrices of each type.
seurat_object= CreateSeuratObject(counts = data$Gene Expression, project = "10X_Seurat")
seurat_object[["CMO"]] = CreateAssayObject(counts = data$Multiplexing Capture)
head([e-mail])
orig.ident nCount_RNA nFeature_RNA nCount_CMO nFeature_CMO
AAACCCAAGAAACACT-1 10X_Seurat 1 1 0 0
AAACCCAAGAAACCAT-1 10X_Seurat 1 1 0 0
AAACCCAAGAAACCCA-1 10X_Seurat 7 7 0 0
AAACCCAAGAAACCCG-1 10X_Seurat 0 0 1 1
AAACCCAAGAAACCTG-1 10X_Seurat 1 1 0 0
AAACCCAAGAAACTAC-1 10X_Seurat 3 3 0 0
library(data.table)
data.table 1.14.2 using 8 threads (see ?getDTthreads). Latest news: r-datatable.com
head(seurat_object@assays)
$RNA
Assay data with 32285 features for 2893103 cells
First 10 features:
Xkr4, Gm1992, Gm19938, Gm37381, Rp1, Sox17, Gm37587, Gm37323, Mrpl15, Lypla1

$CMO
Assay data with 12 features for 2893103 cells
First 10 features:
CMO301, CMO302, CMO303, CMO304, CMO305, CMO306, CMO307, CMO308, CMO309, CMO310

cells <- fread("assignment_confidence_table.csv",select = c("Barcodes"))

cells[1:5]
Barcodes
1: AAACCCAAGACTCAAA-1
2: AAACCCAAGATACATG-1
3: AAACCCAAGATAGCAT-1
4: AAACCCAAGCCTGTCG-1
5: AAACCCAAGTGGCCTC-1
seurat_object_use <- subset(seurat_object, cells = cells$Barcodes)
DefaultAssay(seurat_object_use)<-"CMO"
seurat_object_use <- subset(seurat_object_use, features = c("CMO309","CMO310", "CMO311", "CMO312"))
seurat_object_use = subset(x = seurat_object_use, subset = nCount_CMO > 0)
seurat_object_norm <- NormalizeData(seurat_object_use, assay = "CMO", normalization.method = "CLR")
Normalizing across features
|++++++++++++++++++++++++++++++++++++++++++++++++++| 100% elapsed=00s
seurat_object_demux <- HTODemux(seurat_object_norm, assay = "CMO", positive.quantile = 0.99)
Cutoff for CMO309 : 1285 reads
Cutoff for CMO310 : 966 reads
Cutoff for CMO311 : 1028 reads
Cutoff for CMO312 : 1056 reads
table(seurat_object_demux$CMO_classification.global)

Doublet Negative Singlet
880 1599 9862

Idents(seurat_object_demux) <- "CMO_maxID"

RidgePlot(seurat_object_demux, assay = "CMO", features = rownames(seurat_object_demux[["CMO"]]

RidgePlot(seurat_object_demux, assay = "CMO", features = rownames(seurat_object_demux[["CMO"]])[1:4], ncol = 2)
Picking joint bandwidth of 0.061
Picking joint bandwidth of 0.0577
Picking joint bandwidth of 0.0592
Picking joint bandwidth of 0.064
FeatureScatter(seurat_object_demux, feature1 = "CMO309", feature2 = "CMO310")
Idents(seurat_object_demux) <- "CMO_classification.global"
VlnPlot(seurat_object_demux, features = "nCount_RNA", pt.size = 0.1, log = TRUE)
seurat_object_demux_subset<- subset(seurat_object_demux, idents = "Negative", invert = TRUE)

cmo.dist.mtx <- as.matrix(dist(t(GetAssayData(object = seurat_object_demux_subset, assay = "CMO"))))
seurat_object_demux_subset <- RunTSNE(seurat_object_demux_subset, distance.matrix = cmo.dist.mtx, perplexity = 100)
Warning: Adding a command log without an assay associated with it
DimPlot(seurat_object_demux_subset)
HTOHeatmap(seurat_object_demux, assay = "CMO", ncells = 7000)
Warning message:
guides(<scale> = FALSE) is deprecated. Please use guides(<scale> = "none") instead.
seurat.singlet <- subset(seurat_object_demux, idents = "Singlet")
seurat.singlet <- FindVariableFeatures(seurat.singlet, selection.method = "mean.var.plot")
Calculating gene means
0% 10 20 30 40 50 60 70 80 90 100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Calculating gene variance to mean ratios
0% 10 20 30 40 50 60 70 80 90 100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
**seurat.singlet <- ScaleData(seurat.singlet, features = VariableFeatures(seurat.singlet))
Centering and scaling data matrix
|===========================================================================================================| 100%
seurat.singlet <- RunPCA(seurat.singlet, features = VariableFeatures(seurat.singlet))
Error in irlba(A = t(x = object), nv = npcs, ...) :
max(nu, nv) must be positive
ElbowPlot(seurat.singlet)
Error: Cannot find 'pca' in this Seurat object
ElbowPlot(seurat_object_demux)
Error: Cannot find 'pca' in this Seurat object
View(seurat.singlet)
seurat.singlet <- RunPCA(seurat.singlet, features = VariableFeatures(object = seurat.singlet))
Error in irlba(A = t(x = object), nv = npcs, ...) :
max(nu, nv) must be positive
print(seurat.singlet[["pca"]], dims = 1:5, nfeatures = 5)
Error: Cannot find 'pca' in this Seurat object**

Answer 1

[gorizwango]

Hi shweta, I faced the same problem as you did but I got out of the trouble with the help of a bioinformatician. Apparently, the PCA is absent in your seurat object. You might have missed to run ScaleData, RunPCA and RunUMAP on the integrated data. Following commands may help after you create your integrated object:

seu_int <- Seurat::ScaleData(seu_int)

seu_int <- Seurat::RunPCA(seu_int, npcs = 30)

seu_int <- Seurat::RunUMAP(seu_int, reduction = "pca", dims = 1:25)