STACAS semi-supervised integration completes but during anchor finding throws Error in if (totalCols == 0) return(NULL) : argument is of length zero

[lcoletto]

Hi,
I am running STACAS on a large Seurat object (≈281k cells) with some batches and very unbalanced sample sizes and I would like to clarify whether the behaviour I observe is expected or indicates a problem in the integration.
Below I summarize three different approaches I tried. In all cases, the pipeline runs to completion, but in some cases STACAS prints an error during anchor finding, even though the integration and downstream UMAP are produced.

Dataset / setup
Seurat / SeuratObject recently updated to SeuratV5, RNA assay, log-normalized
Total cells of the object: ~281,600
Number of samples (orig.ident, used as batch): 57
Highly unbalanced batches
11/57 samples have < 1,000 cells
Smallest sample: 265 cells
Largest sample: 13,669 cells
anchor.features = 1000

Cell labels used for semi-supervised mode
"clusters" metadata contains ~38 annotated clusters
Annotated cells: ~38,800
Unannotated cells (NA): 242,828 (majority of the dataset)

Method 1 – Stepwise STACAS, semisupervised, ndim = 28
ndim was chosen based on PCA variance (95% cumulative variance).

library(STACAS)
nfeatures = 1000
ndim = 28
obj.list <- SplitObject(All1, split.by = "orig.ident")
for (n in 1:length(obj.list)) {
  print(n)
  print(obj.list[[n]])
  Idents(obj.list[[n]]) <- "clusters"
}

stacas_anchors <- FindAnchors.STACAS(obj.list, 
                                     anchor.features = nfeatures,
                                     dims = 1:ndim, 
                                     cell.labels = "clusters")

st1 <- SampleTree.STACAS(
  anchorset = stacas_anchors,
  obj.names = names(obj.list))
object_integrated <- IntegrateData.STACAS(stacas_anchors,
                                          sample.tree = st1,
                                          dims=1:ndim)

object_integrated <- object_integrated %>% ScaleData() %>%
  RunPCA(npcs=28) %>% RunUMAP(dims=1:28)

This finishes successfully, but during FindAnchors.STACAS I observe errors like:
Error in if (totalCols == 0) return(NULL) : argument is of length zero
The pipeline does not stop and produces an integrated object.

Method 2: One-liner Run.STACAS, semi-supervised, ndim = 20

library(STACAS)
nfeatures = 1000
ndim = 20
Idents(All1) = "clusters"
object_integrated1 <- All1 %>% SplitObject(split.by = "orig.ident") %>%Run.STACAS(dims = 1:ndim, anchor.features = nfeatures, cell.labels = "clusters") %>% RunUMAP(dims = 1:ndim) 

This also finishes to the end, but I again see the same error message during the run:
Warning: sparse->dense coercion: allocating vector of size 1.0 GiBWarning: pseudoinverse used at -2.2162Warning: neighborhood radius 0.30103Warning: reciprocal condition number 1.4523e-14Preparing PCA embeddings for objects...
|++++++++++++++++++++++++++++++++++++++++++++++++++| 100% elapsed=44s
|+++++++ | 13% ~07h 26m 18s Error in if (totalCols == 0) return(NULL) : argument is of length zero
|++++++++++++++++++++++++++++++++++++++++++++++++++| 100% elapsed=06h 03m 47s

Method 3: One-liner Run.STACAS, unsupervised, ndim = 20

object_integrated2 <- All1 %>% SplitObject(split.by = "orig.ident") %>%Run.STACAS(dims = 1:ndim, anchor.features = nfeatures) %>% RunUMAP(dims = 1:ndim) 

This version:
finishes without explicit errors
produces an integrated object and UMAP

Questions
Is it expected that:
STACAS completes even when FindAnchors.STACAS encounters cases where totalCols == 0?
When using cell.labels, does this error indicate that some batch pairs have no compatible anchors and are effectively skipped?
Is there a recommended way to diagnose which datasets or batch pairs fail to form anchors?
For large, heterogeneous datasets, is semi-supervised STACAS still recommended, or should the unsupervised mode be preferred?

Thank you very much for your help, and for developing STACAS!

[mass-a]

Hi @lcoletto !
thanks for the detailed report.

I am not sure what triggers the error, but I suspect it has to do with some of the smaller samples, for which the method may fail to find integration anchors. To test whether this is the case, I would suggest to re-run the method without the smallest samples (e.g. by dropping the 11 samples with <1000 cells). When you only have 265 cells it may be difficult for the method to find integration anchors, so these small samples may actually contribute more noise than anything...

-> For large, heterogeneous datasets, is semi-supervised STACAS still recommended, or should the unsupervised mode be preferred? Yes, semi-supervised integration is particularly useful for hetereogeneous datasets: the supervised labels can be very helpful in guiding batch-effect correction based on the cell types that are in common between a pair of samples. Always keep in mind, though, that you need at least one cell type in common between a given pair of samples to be able to estimate batch effect correction vectors.

Best
-m

[lcoletto]

Hi @mass-a, thanks a lot for the clarification and the quick reply.

Just to make sure I understand correctly: when you mention that at least one cell type in common between a given pair of samples is required to estimate batch-effect correction vectors, does this mean that, in the semi-supervised setting, each sample needs to share at least one pre-annotated population?

If so, my plan would be to:
assign a coarse, easily identifiable shared population pre-integration (to be used as an anchor),
re-run the integration both including all samples and after excluding the smallest ones (<1000 cells), to test whether the error is driven by limited cell numbers

Thanks again for the helpful suggestions.
Best,
Lavi

[mass-a]

No, I didn't mean that you need a pre-annotated cell type in common between any given pair of samples. Apologies if I have been unclear.

My comment was about having at least a subset of cells of the same type (whether they are annotated or not) between any given pair of samples. For instance, if dataset A only contains macrophages and dataset B only contains B cells, there is no hope of disentangling biological from technical differences and estimating/correcting batch effects. I believe this is true for any scRNA-seq integration algorithm. In semi-supervised STACAS, unannotated cells (marked as NA) are free to form anchors with annotated cells or other NA cells; the semi-supervised labels are simply a way of reducing the space of possible anchors and making life easier for the algorithm. I would recommend to only assign cell type labels you are confident about, and leave the rest to NA.

I re-iterate the suggestion of excluding small samples (e.g. <500 or <1000) to test whether the error stems from the number of cells. You may also find that the integrated space is not very different with or without these samples, so you may decide to simply continue the analysis without the small samples. But let's see what you get :)

Best
-m

[lcoletto]

Hi @mass-a, a quick update.

You were right, re-running the semi-supervised integration after removing the samples with <1000 cells resolved the error, and the pipeline now completes without issues. For the actual analysis I will still need to include those smaller samples, but visually the integrations look very similar when grouping by patient, so this does not seem to introduce obvious artifacts.

Thanks again for the very helpful guidance!
Best,
Lavi