Problem running scPred

Question

[Nitin]

I am trying to predict the cell type using the scPred package. I use Posit R-Studio IDE on the AWS sagemaker.

I am not able to get any output when running the following command

# load the data 
alldata <- readRDS("data/results/run4.rds")

ctrl = alldata[, alldata$orig.ident == 'ctrl']

# set active assay to RNA
[e-mail] = 'RNA' 
ctrl


reference <- scPred::pbmc_1
reference


reference <- reference %>%
  NormalizeData() %>%
  FindVariableFeatures() %>%
  ScaleData() %>%
  RunPCA(verbose = F) %>%
  RunUMAP(dims = 1:30)


#Set the identity resolution 0.5
ctrl <- SetIdent(ctrl, value = "CCA_snn_res.0.5")

ctrl <- ctrl %>%
  NormalizeData() %>%
  FindVariableFeatures() %>%
  ScaleData() %>%
  RunPCA(verbose = F) %>%
  RunUMAP(dims = 1:30)

# Transfer Label 
transfer.anchors <- FindTransferAnchors(reference = reference, query = ctrl, 
    dims = 1:30)
predictions <- TransferData(anchorset = transfer.anchors, refdata = reference$cell_type, 
    dims = 1:30)
ctrl <- AddMetaData(object = ctrl, metadata = predictions)

## scPred 
reference <- getFeatureSpace(reference, "cell_type")
reference <- trainModel(reference)

get_scpred(reference)

Till here everything works fine, but the following code fails

ctrl <- scPredict(ctrl, reference)

OUTPUT

Matching reference with new dataset...
     2000 features present in reference loadings
     1773 features shared between reference and new dataset
     88.65% of features in the reference are present in new dataset
Error in `GetAssayData()`:
! `assay` must be one of "RNA", not "data".
Backtrace:
 1. scPred::scPredict(ctrl, reference)
 2. scPred::project_query(...)
 4. SeuratObject:::GetAssayData.Seurat(assay = "data")
 Error in GetAssayData(new, "data") :

Please help

[mazegriff]

Hello,

I had the exact same issue and found there is a bug in the package. Specifically, in scPred/R/project_query.R on line 87. The line has the following code:

new_data <- GetAssayData(new, "data")[shared_features,]

The Seurat function GetAssayData() requires the third positional argument to be stated as slot="data"; otherwise, its second positional argument assay expects a suitable value like "RNA" or another type of assay (and will throw the error you cited). Here is the correction of line 87 in the script:

new_data <- GetAssayData(new, slot="data")[shared_features,]

When I cloned the repo to the remote machine I perform work on, I customized the package by modifying this line and successfully implemented it in my R runtime environment (R 4.3).

I've notified the maintainer of the package of this bug.

Hope this helps.

Maze

[tfguinan]

This is correct, after updating Seurat to v5 the error occurs due to GetAssayData() now requiring a named parameter (and the second position now defaulting to 'assay').

Worth noting the 'slot' param is now deprecated as per the SeuratObject documentation. So the pull request uses layer="data" for continued functionality.

scPred with the updated call for the v5 function can be installed by referencing the specific pull request:

devtools::install_github(repo="powellgenomicslab/scPred", ref="9f407b7436f40d44224a5976a94cc6815c6e837f", force = TRUE)

[Ram]

[e-mail] = 'RNA'

Are you sure this is the right way to set the active assay?

[Nitin]

yes, also I tried with

DefaultAssay(object = ctrl) <- "RNA"

[Ram]

Can you run these commands before setting the default assay and show us the output:

Assays(ctrl)
DefaultAssay(ctrl)
DefaultAssay(ctrl) <- 'RNA'
DefaultAssay(ctrl)

[Nitin]

Here is the output

Assays(ctrl)

An object of class "SimpleAssays"
Slot "data":
List of length 1


DefaultAssay(ctrl)
[1] "CCA"

DefaultAssay(ctrl) <- 'RNA'

DefaultAssay(ctrl)
[1] "RNA"

[Ram]

That does not make sense to me. You don't have an RNA assay, yet you're assigning it as the default assay. Does that add up to you?

[Nitin]

yes, but when check the

ctrl

An object of class Seurat 
38472 features across 29125 samples within 2 assays 
Active assay: RNA (36472 features, 2000 variable features)
 3 layers present: counts, data, scale.data
 1 other assay present: CCA
 7 dimensional reductions calculated: pca, umap, tsne, harmony, umap_harmony, scanorama, umap_scanorama

[Ram]

Oh, OK. I guess I was wrong. What does Assays(ctrl) give you? Also, what does class(ctrl[["RNA"]]) say?

[Nitin]

NOTE - Also I tried to run the scPRED as per the link https://powellgenomicslab.github.io/scPred/articles/introduction.html I have the same issue at the same step.

Here are the outputs you asked for

Assays(ctrl)

An object of class "SimpleAssays"
Slot "data":
List of length 1


class(ctrl[["RNA"]])

[1] "Assay"
attr(,"package")
[1] "SeuratObject"

[Ram]

You're going to need to consult a bioinformatician close to you or wait for a Seurat expert here to help you out.

Answer 1

[vtm]

I did solve the problem using these steps

1. Run the code

trace(scPred::project_query, edit=TRUE)

2. Modify the code as below in the popped up window and save

new_data <- GetAssayData(new, layer = "data")[shared_features, ]