Transfer harmony-integrated scRNA-seq data to scanpy

Question

[fifty_fifty]

I have integrated two scRNA-seq datasets using harmony method. Unlike Seurat integration, Harmony just adds extra embeddings for further clustering and other analysis. However, the data (raw counts, normalized, and scaled slots) is the same as in unintegrated data. So, when I transfer the normalized counts to scanpy for other analysis (ikarus) obviously the data behaves as unintegrated. Is there any way to recompute normalized counts using harmony embeddings so the data becomes transferable to scanpy as an integrated dataset?

[Madiha]

Hi, I am not sure if you still remmeber the solution you followed that time. But I am stuck with similar problem. scVI is not working good fro my data and if i integrate and batch correct using harmony it shows good integration on umap in R but when I transfer data to python adata umap shows batch effects and unintegrated like data. Wondering what you did to export Harmony integration properly.

Thanks

Answer 1

[Radu Tanasa]

As you well noticed, harmony creates an embedding, it does not generate integrated normalized counts. If you need that, I suggest you use scVI, which does that very well. In integration benchmarks, scVI generally performed better. Here's a simple example. Make sure you have the "batch" observation, which marks the cells from each dataset.

sc.pp.highly_variable_genes(
  adata,
  flavor="seurat_v3",
  n_top_genes=2000,
  layer="counts",
  batch_key="batch"
)
adata.layers['counts'] = adata.X.copy()
scvi.model.SCVI.setup_anndata(adata, layer="counts", batch_key="batch")
model = scvi.model.SCVI(adata, n_layers=2, n_latent=30, gene_likelihood="nb")
adata.obsm['X_scVI'] = model.get_latent_representation()
adata.layers['scvi_normalized'] = model.get_normalized_expression(library_size = 1e4)

Hopefully, you can then use the scVI normalized expression for further downstream analysis.