Regress Out using all genes

[Khalid-Usman]

Why can't I use regress_out function for scRNA-seq data without applying highly_variable_genes. Also I think regress_out function should be before highly_variable_genes, because in this way we can first remove batch effect and then select important genes.

sc.pp.highly_variable_genes(adata, min_mean=0.0125, max_mean=3, min_disp=0.5)
sc.pl.highly_variable_genes(adata)
#adata = adata[:, adata.var['highly_variable']]
print(adata)

Regressing

sc.pp.regress_out(adata, ['n_counts', 'percent_mito'])

[LuckyMD]

Hi @Khalid-Usman,

Regressing out should indeed be performed before highly variable gene selection. This was not in the original scRNA-seq tutorials from Seurat and Scanpy though. If you're interested in a current best-practices tutorial (based on scanpy, but also including R tools), you can find it here.

The reason it might not have been done on all genes initially is for speed.

[GuidoBarzaghi]

Hi everyone!

I was wondering about the same issue. Would you then suggest to regress after subsetting HVGs (for speed reasons) and then re-searching and re-subsetting HVGs after the unwanted source of variation is corrected for? In such a way, one would inevitably loose some interesting genes (due to the first HVGs subsetting), but a cleaner signal would be obtained at the end.

[LuckyMD]

I don't think I would do HVG selection twice. I reckon it's often not such a time-taking process to regress out across all genes. Just do it once, and then do HVG selection afterwards. However, with the current setup this is difficult as sc.pp.regress_out also removes the offset, and thus has 0 mean per gene. This prevents sc.pp.highly_variable from binning the genes by mean expression (see #722).

[Hrovatin]

Is HVG function designed to be done after regress out?

# Regress out process of interest from expression data
sc.pp.regress_out(adata_b_rn_sub2, keys='LogReg_decision') 
# Find HVGs (across samples, not per sample as samples are very different in terms of beta bio)
sc.pp.highly_variable_genes(adata_b_rn_sub2, flavor='cell_ranger', 
                            batch_key=None, n_top_genes =2000)
print('\n','Number of highly variable genes: {:d}'.format(
    np.sum(adata_b_rn_sub2.var['highly_variable'])))
rcParams['figure.figsize']=(10,5)
sc.pl.highly_variable_genes(adata_b_rn_sub2)

If I do scale before HVG

[LuckyMD]

As I wrote above, our standard hvg selection function bins genes by mean, and thus doesn't work well with the current implementation of regress_out. You would have to add the offset again to do this properly, and not just use the residuals.

[Hrovatin]

Oh I somehow overlooked the last sentence.

[VivianeSchulz]

I have the same problem and I do not know how to 'add the offset again'. Could you maybe tell me how this is done?

Just for context: I am regressing out cell cycle and mitochondrial genes, after which I scale and do PCA. I then analyse my dataset, but when I want to sub-cluster, I want to run hvg selection again. Since that is not possible with the current implementation of regress_out, I am trying to find a solution for this.