Integrating multiple objects with IntegrateEmbeddings() erases all previous reductions

[hannanw]

I have four different 10x multiome seurat object (one for each experiment). I am able to integrate the RNA assay without much issue.

merge_obj <- merge(expt1, c(expt2,expt3,expt4))
DefaultAssay(merge_obj) <- "SCT"
obj.list <- SplitObject(merge_obj, split.by = "orig.ident")
obj.features <- SelectIntegrationFeatures(obj.list, nfeatures = 3000)
VariableFeatures(merge_obj[["SCT"]]) <- obj.features
merge_obj <- RunPCA(merge_obj, dims=1:50, verbose = FALSE)
merge_obj_integration <- IntegrateLayers(object = merge_obj, method = CCAIntegration, normalization.method = "SCT", verbose = F)
merge_obj_integration <- RunPCA(merge_obj_integration, dims=1:50, verbose = FALSE)
merge_obj_integration <- FindNeighbors(merge_obj_integration, reduction = "integrated.dr", dims = 1:16)
merge_obj_integration <- FindClusters(merge_obj_integration, res=0.2,method='igraph',algorithm=4)

>merge_obj_integration
An object of class Seurat
243934 features across 23666 samples within 3 assays
Active assay: SCT (23167 features, 2960 variable features)
 3 layers present: counts, data, scale.data
 2 other assays present: RNA, ATAC
 2 dimensional reductions calculated: pca, integrated.dr

However, when I try to do the same integration for the ATAC assay IntegrateEmbeddings() removes all prior reductions. I followed the step from the Signac scATAC-seq data integration here.

DefaultAssay(merge_obj_integration) <- "ATAC"
merge_obj_integration <- FindTopFeatures(merge_obj_integration, min.cutoff = 10)
merge_obj_integration <- RunTFIDF(merge_obj_integration)
merge_obj_integration <- RunSVD(merge_obj_integration)
merge_obj_integration <- RunUMAP(merge_obj_integration, reduction = "lsi", dims = 2:50)

> merge_obj_integration
An object of class Seurat
243934 features across 23666 samples within 3 assays
Active assay: ATAC (184166 features, 183851 variable features)
 2 layers present: counts, data
 2 other assays present: RNA, SCT
 5 dimensional reductions calculated: pca, integrated.dr, integrated_RNA_umap, lsi, umap

atac.obj.list <- SplitObject(merge_obj_integration, split.by = "orig.ident")

integration.anchors <- FindIntegrationAnchors(
  object.list = atac.obj.list,
  anchor.features = rownames(merge_obj_integration),
  reduction = "rlsi",
  dims = 2:50
)

merge_obj_integration_ATAC <- IntegrateEmbeddings(
  anchorset = integration.anchors,
  reductions = merge_obj_integration[["lsi"]],
  k.weight = 50,
  new.reduction.name = "integrated_lsi",
  dims.to.integrate = 2:50
)

> merge_obj_integration_ATAC
An object of class Seurat
243934 features across 23666 samples within 3 assays
Active assay: ATAC (184166 features, 0 variable features)
 2 layers present: counts, data
 2 other assays present: RNA, SCT
 1 dimensional reduction calculated: integrated_lsi

Am I missing or doing something wrong? Any help/advice would be greatly appreciated!

> sessionInfo()
R version 4.3.2 (2023-10-31)
Platform: x86_64-conda-linux-gnu (64-bit)
Running under: Ubuntu 18.04.6 LTS

attached base packages:
[1] parallel  stats4    stats     graphics  grDevices utils     datasets
[8] methods   base

other attached packages:
 [1] ROCR_1.0-11                       KernSmooth_2.23-22
 [3] sctransform_0.4.1                 fields_15.2
 [5] viridisLite_0.4.2                 spam_2.10-0
 [7] leiden_0.4.3.1                    reticulate_1.35.0
 [9] scDblFinder_1.16.0                SingleCellExperiment_1.24.0
[11] SummarizedExperiment_1.32.0       MatrixGenerics_1.14.0
[13] matrixStats_1.2.0                 DoubletFinder_2.0.4
[15] glmGamPoi_1.14.3                  hdf5r_1.3.9
[17] BSgenome.Hsapiens.UCSC.hg38_1.4.5 BSgenome_1.70.2
[19] rtracklayer_1.62.0                BiocIO_1.12.0
[21] Biostrings_2.70.2                 XVector_0.42.0
[23] EnsDb.Hsapiens.v86_2.99.0         ensembldb_2.26.0
[25] AnnotationFilter_1.26.0           GenomicFeatures_1.54.3
[27] AnnotationDbi_1.64.1              Biobase_2.62.0
[29] GenomicRanges_1.54.1              GenomeInfoDb_1.38.6
[31] IRanges_2.36.0                    S4Vectors_0.40.2
[33] BiocGenerics_0.48.1               ggforce_0.4.2
[35] qlcMatrix_0.9.7                   sparsesvd_0.2-2
[37] slam_0.1-50                       Matrix_1.6-5
[39] Signac_1.12.0                     ggplot2_3.4.4
[41] dplyr_1.1.4                       SeuratWrappers_0.3.4
[43] Seurat_5.0.1                      SeuratObject_5.0.1
[45] sp_2.1-3

[timoast]

Here you're actually looking at a different object than the one you computed the other reductions for. When you split and then merge an object (as is done in the integration function), the reductions are removed. The easiest way around this would be to just add the integrated LSI reduction to the original merged object merge_obj_integration