satijalab / seurat

R toolkit for single cell genomics
http://www.satijalab.org/seurat
Other
2.26k stars 907 forks source link

Error: Given reference assay (SCT) has 5 reference sct models. Please provide a reference assay with a single reference sct model. #8566

Open xmliu01 opened 7 months ago

xmliu01 commented 7 months ago

I have two data sets of single-cell transcriptomes: snrna and sprna, where snrna was obtained under Seurat v4.3 and sprna was obtained under Seurat v5. I use the same code, as shown below.

W3.25_1.count=Read10X("./data/W3.25_1/",gene.column=1)
W3.25_2.count=Read10X("./data/W3.25_2/",gene.column=1)
W3.25_3.count=Read10X("./data/W3.25_3/",gene.column=1)
W4_1.count=Read10X("./data/W4_1/",gene.column=1)
W4_2.count=Read10X("./data/W4_2/",gene.column=1)

W3.25_1=CreateSeuratObject(counts=W3.25_1.count,project = "W3.25_1",min.cells=3,min.features=100)
W3.25_2=CreateSeuratObject(counts=W3.25_2.count,project = "W3.25_2",min.cells=3,min.features=100)
W3.25_3=CreateSeuratObject(counts=W3.25_3.count,project = "W3.25_3",min.cells=3,min.features=100)
W4_1=CreateSeuratObject(counts=W4_1.count,project = "W4_1",min.cells=3,min.features=100)
W4_2=CreateSeuratObject(counts=W4_2.count,project = "W4_2",min.cells=3,min.features=100)

all=merge(W3.25_1,list(W3.25_2,W3.25_3,W4_1,W4_2),add.cell.ids=c("W3.25_1","W3.25_2","W3.25_3","W4_1","W4_2"),project="all")
all=subset(all,nFeature_RNA > 100 & nFeature_RNA<10000 & nCount_RNA>100 & nFeature_RNA<10000)

all=SCTransform(all,variable.features.n = 2000)
all=RunPCA(all,features = VariableFeatures(object = all))
all=RunUMAP(all,dims=1:10,return.model=T)
all=FindNeighbors(all,dims=1:10)
all=FindClusters(all,resolution = 0.5)

However, only one SCT model was obtained in Seurat v4.3, while the sprna obtained in Seurat v5 had five SCT models. When I try to compare these two datasets, I get the following error with FindTransferAnchors.

Error: Given reference assay (SCT) has 5 reference sct models. Please provide a reference assay with a  single reference sct model.

How do I get one SCT model in Seurat v5?

ollieeknight commented 7 months ago

with seurat 5, you need to run JoinLayers first, like

all <- JoinLayers(all, assay = 'SCT')
xmliu01 commented 7 months ago

with seurat 5, you need to run JoinLayers first, like

all <- JoinLayers(all, assay = 'SCT')

Thanks, I tried it, but it not work.

ollieeknight commented 7 months ago

is your assay named SCT, or RNA? you can also just try JoinLayers(all)

saketkc commented 7 months ago

@xmliu01 Can you paste your full code here?

JHCCoder commented 6 months ago

Is FindTransferAnchors available supported for SCTransformed, integrated data? Here is a minimally reproducible problem with sample data.

library(Seurat)
library(SeuratData)
library(ggplot2)
panc8 <- LoadData("panc8")
table(panc8$tech)

Reference processing

pancreas.ref <- subset(panc8, tech %in% c("celseq2", "smartseq2"))
pancreas.ref[["RNA"]] <- split(pancreas.ref[["RNA"]], f = pancreas.ref$tech)
## Sctransform normalize
pancreas.ref <- SCTransform(pancreas.ref)
pancreas.ref <- RunPCA(pancreas.ref)
pancreas.ref <- RunUMAP(pancreas.ref, dims = 1:30)

## Plot 
options(repr.plot.height = 8, repr.plot.width = 17)
DimPlot(pancreas.ref, reduction = "umap", group.by = c("tech", "celltype"))

download - 2024-03-10T000707 019

# integrate datasets
pancreas.ref <- IntegrateLayers(object = pancreas.ref, method = CCAIntegration, normalization.method = "SCT", verbose = F)

This is what pancreas.ref looks like: An object of class Seurat 55515 features across 4679 samples within 2 assays Active assay: SCT (21152 features, 3000 variable features) 3 layers present: counts, data, scale.data 1 other assay present: RNA 2 dimensional reductions calculated: pca, umap

## Visualize the effect of integration on embedding
pancreas.ref <- RunUMAP(pancreas.ref, dims = 1:30, reduction = "integrated.dr")
DimPlot(pancreas.ref, reduction = "umap", group.by = c("tech", "celltype"))

download - 2024-03-10T001132 428 Looks effectively batch corrected, great.

This is what pancreas.ref looks like:

An object of class Seurat 
55515 features across 4679 samples within 2 assays 
Active assay: SCT (21152 features, 3000 variable features)
 3 layers present: counts, data, scale.data
 1 other assay present: RNA
 3 dimensional reductions calculated: pca, umap, integrated.dr

And for query processing same thing:

# select two technologies for the query datasets
pancreas.query <- subset(panc8, tech %in% c("fluidigmc1", "celseq"))
# we will use data from 2 technologies for the reference
pancreas.query[["RNA"]] <- split(pancreas.query[["RNA"]], f = pancreas.query$tech)
pancreas.query <- SCTransform(pancreas.query)

This is what query data looks like at this point

An object of class Seurat 
55828 features across 1642 samples within 2 assays 
Active assay: SCT (21465 features, 3000 variable features)
 3 layers present: counts, data, scale.data
 1 other assay present: RNA

Then for the part of FindTransferAnchor here are the following code:

pancreas.anchors <- FindTransferAnchors(reference = pancreas.ref,
                                        query = pancreas.query, 
                                        normalization.method="SCT",
                                        reference.assay="SCT",
                                        dims = 1:30,
                                        reference.reduction = "pca")

And the error:

Error: Given reference assay (SCT) has 2 reference sct models. Please provide a reference assay with a  single reference sct model.
Traceback:

1. FindTransferAnchors(reference = pancreas.ref, query = pancreas.query, 
 .     normalization.method = "SCT", reference.assay = "SCT", dims = 1:30, 
 .     reference.reduction = "pca")
2. ValidateParams_FindTransferAnchors(reference = reference, query = query, 
 .     normalization.method = normalization.method, recompute.residuals = recompute.residuals, 
 .     reference.assay = reference.assay, reference.neighbors = reference.neighbors, 
 .     query.assay = query.assay, reduction = reduction, reference.reduction = reference.reduction, 
 .     project.query = project.query, features = features, scale = scale, 
 .     npcs = npcs, l2.norm = l2.norm, dims = dims, k.anchor = k.anchor, 
 .     k.filter = k.filter, k.score = k.score, max.features = max.features, 
 .     nn.method = nn.method, n.trees = n.trees, eps = eps, approx.pca = approx.pca, 
 .     mapping.score.k = mapping.score.k, verbose = verbose)
3. stop("Given reference assay (", reference.assay, ") has ", reference.model.num, 
 .     " reference sct models. Please provide a reference assay with a ", 
 .     " single reference sct model.", call. = FALSE)

There is the two reference SCTmodels? The SCT assay seems to only indicate one model? If you can identify the issue here, any help would be appreciated...

JABioinf commented 6 months ago

I have the same problem since I updated Seurat to 5.0.2 Notably there seems to have been a major change with either SCTransform or IntegrateLayers() since the later doesn't generate an "integrated" assay anymore and only a dimensionality reduction "integrated.dr". How does that impact FindTransferAnchors when applying SCT? Similarly since the update, even FindTransferAnchors codes that previously worked using the RNA normalization pipeline now run indefinitely without throwing error. Likewise IntegrateLayers did not generated "integrated" assay in that case. Was that expected from those updates? Could the Seurat team update or verify how the integration process should be done both for SCT and RNA lognorm approach?

igrabski commented 3 months ago

Thanks for providing all this detail -- we are currently investigating this problem and hope to have a solution soon.