# Preparing data for totalVI

**URL:** <https://discourse.scverse.org/t/preparing-data-for-totalvi/149>\
**Category:** scvi-tools\
**Tags:** totalvi, preprocessing\
**Created:** [July 28, 2021, 11:19am UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149 "2021-07-28T11:19:01Z")\
**Posts on this page:** 12\
**Page:** 1

<div class="post-metadata">

**Author:** ![michael101](https://avatars.discourse-cdn.com/v4/letter/m/df705f/32.png) [@michael101](https://discourse.scverse.org/u/michael101)\
**Post date:** [July 28, 2021, 11:19am UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149/1 "2021-07-28T11:19:01Z")

</div>

Hello,

I am hoping to ask about TOTALVI. I have followed the steps in the tutorial, however when I come to do vae.train() I get the following error which I do not understand:

```auto
GPU available: True, used: True
TPU available: False, using: 0 TPU cores
IPU available: False, using: 0 IPUs
/usr/local/lib/python3.7/dist-packages/pytorch_lightning/core/datamodule.py:424: LightningDeprecationWarning: DataModule.setup has already been called, so it will not be called again. In v1.6 this behavior will change to always call DataModule.setup.
  f"DataModule.{name} has already been called, so it will not be called again. "
LOCAL_RANK: 0 - CUDA_VISIBLE_DEVICES: [0]
/usr/local/lib/python3.7/dist-packages/pytorch_lightning/trainer/data_loading.py:323: UserWarning: The number of training samples (34) is smaller than the logging interval Trainer(log_every_n_steps=50). Set a lower value for log_every_n_steps if you want to see logs for the training epoch.
  f"The number of training samples ({self.num_training_batches}) is smaller than the logging interval"
Epoch 1/400: -0%| | -1/400 [00:00<?, ?it/s]/usr/local/lib/python3.7/dist-packages/tqdm/std.py:538: TqdmWarning: clamping frac to range [0, 1]
  colour=colour)
Epoch 1/400: -0%| | -1/400 [00:00<?, ?it/s]
---------------------------------------------------------------------------
KeyError Traceback (most recent call last)
<ipython-input-70-511a56042b55> in <module>()
----> 1 vae.train()

15 frames
/usr/local/lib/python3.7/dist-packages/scvi/train/_trainingplans.py in training_epoch_end(self, outputs)
    133 n_obs, elbo, rec_loss, kl_local = 0, 0, 0, 0
    134 for tensors in outputs:
--> 135 elbo += tensors["reconstruction_loss_sum"] + tensors["kl_local_sum"]
    136 rec_loss += tensors["reconstruction_loss_sum"]
    137 kl_local += tensors["kl_local_sum"]

KeyError: 'reconstruction_loss_sum'

```

I’ve pretty much followed the instructions to the letter so I’m not really sure why I get this error, or indeed what it means. I wonder if you could please advise? Many thanks in advance.

---

<div class="post-metadata">

**Author:** ![michael101](https://avatars.discourse-cdn.com/v4/letter/m/df705f/32.png) [@michael101](https://discourse.scverse.org/u/michael101)\
**Post date:** [July 28, 2021, 11:20am UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149/2 "2021-07-28T11:20:35Z")

</div>

I should probably add that this is what I get when I do the scvi.data.setup\_anndata

```auto
INFO Using batches from adata.obs["batch"]                                               
INFO No label_key inputted, assuming all cells have same label                           
INFO Using data from adata.layers["counts"]                                              
INFO Computing library size prior per batch                                              
INFO Using protein expression from adata.obsm['protein_expression']                      
INFO Using protein names from columns of adata.obsm['protein_expression']                
INFO Found batches with missing protein expression                                       
INFO Successfully registered anndata object containing 9564 cells, 4000 vars, 2 batches, 
         1 labels, and 137 proteins. Also registered 0 extra categorical covariates and 0    
         extra continuous covariates.                                                        
INFO Please do not further modify adata until model is trained.

```

---

<div class="post-metadata">

**Author:** ![michael101](https://avatars.discourse-cdn.com/v4/letter/m/df705f/32.png) [@michael101](https://discourse.scverse.org/u/michael101)\
**Post date:** [July 28, 2021, 3:40pm UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149/3 "2021-07-28T15:40:06Z")

</div>

UPDATE…

Doing scvi.data.organize\_cite\_seq\_10x(adata) AFTER concatenating the two objects I have seems to then make it work.

Even though in the tutorial this is applied BEFORE concatenation.

Doing it before leaves a full repertoire of 137 proteins. Doing it after reduces this to 95.

If anyone knows why this is that information would be much appreciated.

---

<div class="post-metadata">

**Author:** ![adamgayoso](https://yyz1.discourse-cdn.com/flex035/user_avatar/discourse.scverse.org/adamgayoso/32/100_2.png) [@adamgayoso](https://discourse.scverse.org/u/adamgayoso)\
**Post date:** [July 28, 2021, 5:12pm UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149/4 "2021-07-28T17:12:25Z")

</div>

> [@michael101](#):
>
> Doing scvi.data.organize\_cite\_seq\_10x(adata) AFTER concatenating the two objects I have seems to then make it work.

It’s difficult for me to understand what steps are being taken for data organization and preprocessing, but this can make sense to me.

---

<div class="post-metadata">

**Author:** ![michael101](https://avatars.discourse-cdn.com/v4/letter/m/df705f/32.png) [@michael101](https://discourse.scverse.org/u/michael101)\
**Post date:** [July 28, 2021, 8:10pm UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149/5 "2021-07-28T20:10:47Z")

</div>

Thanks for the reply. The steps undertake so far are:

1. Importing the data
2. Scanpy preprocessing/QC metrics
3. Scrublet

These steps were performed separately for two samples/objects. I then did:

1. Concatenating the two objects
2. Filtering min genes/cells/counts/mito
3. Scanpy tools pca/nearest neighbour/umap
4. Then the rest is as per the TOTALVI tutorial (normalise, log, adata.raw, highly variable genes, scvi.data.organize, scvi.data.setup\_anndata)

Not not really clear why the protein count drops off with this method though compared to doing scvi.data.organize on the two samples separately before concatenating them. Surely we are not losing cells because the QC steps are otherwise the same?

---

<div class="post-metadata">

**Author:** ![adamgayoso](https://yyz1.discourse-cdn.com/flex035/user_avatar/discourse.scverse.org/adamgayoso/32/100_2.png) [@adamgayoso](https://discourse.scverse.org/u/adamgayoso)\
**Post date:** [July 29, 2021, 4:54am UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149/6 "2021-07-29T04:54:35Z")

</div>

A script showing your exact steps would be helpful to diagnose the issue. For example, Scrublet should probably not be run on an anndata object containing both RNA and protein in `adata.X` which is how it’s loaded from a 10x output.

---

<div class="post-metadata">

**Author:** ![michael101](https://avatars.discourse-cdn.com/v4/letter/m/df705f/32.png) [@michael101](https://discourse.scverse.org/u/michael101)\
**Post date:** [July 29, 2021, 3:27pm UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149/7 "2021-07-29T15:27:00Z")

</div>

Thanks, Adam. That may well be part of the problem. Apologies, I had not realised Scrublet should not be run on adata with RNA and protein.

It looks like I cannot upload an ipynb file here, but this is the script I’m working from:

# Install required packages

!pip install scvi-tools -U  
!pip install scanpy  
!pip install scrublet  
!pip install leidenalg  
!pip install scrublet

import sys

#if branch is stable, will install via pypi, else will install from source  
branch = “stable”  
IN\_COLAB = “google.colab” in sys.modules

if IN\_COLAB and branch == “stable”:  
!pip install --quiet scvi-tools[tutorials]  
elif IN\_COLAB and branch != “stable”:  
!pip install --quiet --upgrade jsonschema  
!pip install --quiet git+https://github.com/yoseflab/scvi-tools@$branch#egg=scvi-tools[tutorials]

# Import packages

import scanpy as sc  
import scvi  
import matplotlib.pyplot as plt  
import scrublet as scr

#Scanpy set up  
sc.settings.verbosity = 3 # verbosity: errors (0), warnings (1), info (2), hints (3)

#Mount googledrive  
from google.colab import drive  
drive.mount(‘/content/gdrive/’, force\_remount=True)

#Set working directory  
import os  
os.chdir(‘/content/gdrive/MyDrive/Adenoid\_pilot’)

## Importing data from Ad20

ad20 = sc.read\_10x\_h5(‘/content/gdrive/MyDrive/Adenoid\_pilot/filtered\_feature\_bc\_ad20\_matrix.h5’,  
genome=None,  
gex\_only=False,  
backup\_url=None)

#Organise and check the data  
ad20.var\_names\_make\_unique()  
ad20.raw=ad20  
ad20.var[“mito”] = ad20.var\_names.str.startswith(“MT-”)  
sc.pp.calculate\_qc\_metrics(ad20, qc\_vars=[“mito”], inplace=True)

#scrublet  
scrub = scr.Scrublet(ad20.raw.X)  
ad20.obs[‘doublet\_scores’], ad20.obs[‘predicted\_doublets’] = scrub.scrub\_doublets()  
scrub.plot\_histogram()

sum(ad20.obs[‘predicted\_doublets’])

ad20.obs[‘doublet\_info’] = ad20.obs[“predicted\_doublets”].astype(str)

sc.pl.violin(ad20, ‘n\_genes\_by\_counts’,  
jitter=0.4, groupby = ‘doublet\_info’, rotation=45)

ad20.var[“sample”] = “20”

## Importing data from Ad19

ad19 = sc.read\_10x\_h5(‘/content/gdrive/MyDrive/Adenoid\_pilot/filtered\_feature\_bc\_ad19\_matrix.h5’,  
genome=None,  
gex\_only=False,  
backup\_url=None) # write a cache file for faster subsequent reading

#Organise and check the data  
ad19.var\_names\_make\_unique()  
ad19.raw=ad19  
ad19.var[“mito”] = ad19.var\_names.str.startswith(“MT-”)  
sc.pp.calculate\_qc\_metrics(ad19, qc\_vars=[“mito”], inplace=True)

#scrublet  
scrub = scr.Scrublet(ad19.raw.X)  
ad19.obs[‘doublet\_scores’], ad19.obs[‘predicted\_doublets’] = scrub.scrub\_doublets()  
scrub.plot\_histogram()

sum(ad19.obs[‘predicted\_doublets’])

ad19.obs[‘doublet\_info’] = ad19.obs[“predicted\_doublets”].astype(str)

sc.pl.violin(ad19, ‘n\_genes\_by\_counts’,  
jitter=0.4, groupby = ‘doublet\_info’, rotation=45)

ad19.var[“sample”] = “19”

### Combine objects

adata = ad19.concatenate(ad20)  
adata.raw=adata  
adata.layers[“counts”] = adata.X.copy()

# Filter for min cells/genes/counts

sc.pp.filter\_cells(adata, min\_genes=200)  
sc.pp.filter\_genes(adata, min\_cells=3)  
sc.pp.filter\_genes(adata, min\_counts=1)

#Examine mito genes  
sc.pp.calculate\_qc\_metrics(adata, qc\_vars=[“mito”], inplace=True)

#Plot QC metrics  
sc.pl.violin(adata, [‘n\_genes\_by\_counts’, ‘total\_counts’, ‘pct\_counts\_mito’],  
jitter=0.4, multi\_panel=True)

#Filter mito genes  
adata= adata[adata.obs.pct\_counts\_mito \< 10, :]

sc.tl.pca(adata, svd\_solver=‘arpack’)  
sc.pl.pca\_variance\_ratio(adata, log=True, n\_pcs = 50)  
sc.pp.neighbors(adata, n\_pcs = 50, n\_neighbors = 10)  
sc.tl.umap(adata)

sc.pl.umap(adata, color=“doublet\_info”)

adata=adata[adata.obs[‘doublet\_info’] == “False”]

### Prepare for SCVI modelling

adata.layers[“counts”] = adata.X.copy()  
sc.pp.normalize\_total(adata, target\_sum=1e4)  
sc.pp.log1p(adata)  
adata.raw = adata

sc.pp.highly\_variable\_genes(  
adata,  
n\_top\_genes=4000,  
flavor=“seurat\_v3”,  
batch\_key=“batch”,  
subset=True,  
layer=“counts”  
)

scvi.data.organize\_cite\_seq\_10x(adata) #This moves the protein bit into obsm

scvi.data.setup\_anndata(  
adata,  
layer=“counts”,  
batch\_key=“batch”,  
protein\_expression\_obsm\_key=“protein\_expression”  
)

scvi.data.view\_anndata\_setup(adata)

#Training the model  
vae = scvi.model.TOTALVI(adata, latent\_distribution=“normal”)

vae.train()

plt.plot(vae.history[“elbo\_train”], label=“train”)  
plt.plot(vae.history[“elbo\_validation”], label=“validation”)  
plt.title(“Negative ELBO over training epochs”)  
plt.ylim(0, 1400)  
plt.legend()

---

<div class="post-metadata">

**Author:** ![adamgayoso](https://yyz1.discourse-cdn.com/flex035/user_avatar/discourse.scverse.org/adamgayoso/32/100_2.png) [@adamgayoso](https://discourse.scverse.org/u/adamgayoso)\
**Post date:** [July 29, 2021, 5:48pm UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149/8 "2021-07-29T17:48:50Z")

</div>

Looking into this more, it could be a problem with one of our dependencies having a recent update. I’ll get back to you soon.

---

<div class="post-metadata">

**Author:** ![adamgayoso](https://yyz1.discourse-cdn.com/flex035/user_avatar/discourse.scverse.org/adamgayoso/32/100_2.png) [@adamgayoso](https://discourse.scverse.org/u/adamgayoso)\
**Post date:** [July 29, 2021, 7:28pm UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149/9 "2021-07-29T19:28:06Z")

</div>

We just pushed a new version of scvi-tools that fixes part of these issues.

The concatenation before or after the `organize` function shouldn’t really matter… but I use [concat](https://anndata.readthedocs.io/en/latest/anndata.concat.html) instead of concatenate.

---

<div class="post-metadata">

**Author:** ![michael101](https://avatars.discourse-cdn.com/v4/letter/m/df705f/32.png) [@michael101](https://discourse.scverse.org/u/michael101)\
**Post date:** [July 30, 2021, 11:40am UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149/10 "2021-07-30T11:40:00Z")

</div>

Thanks, Adam. I’ve just tried to re-run it but we’re still seeing a fair bit of protein drop out (137-\>95). I’m not really sure why that is happening. Presumably the hope was that the new version would correct that?

---

<div class="post-metadata">

**Author:** ![adamgayoso](https://yyz1.discourse-cdn.com/flex035/user_avatar/discourse.scverse.org/adamgayoso/32/100_2.png) [@adamgayoso](https://discourse.scverse.org/u/adamgayoso)\
**Post date:** [July 30, 2021, 3:09pm UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149/11 "2021-07-30T15:09:41Z")

</div>

The update fixes this:

```auto
/usr/local/lib/python3.7/dist-packages/scvi/train/_trainingplans.py in training_epoch_end(self, outputs)
    133 n_obs, elbo, rec_loss, kl_local = 0, 0, 0, 0
    134 for tensors in outputs:
--> 135 elbo += tensors["reconstruction_loss_sum"] + tensors["kl_local_sum"]
    136 rec_loss += tensors["reconstruction_loss_sum"]
    137 kl_local += tensors["kl_local_sum"]

```

There could be several other reasons the number of proteins is dropping and it’s probably related to the data manipulation/filtering. This could include the scanpy filtering steps (as proteins are in `adata.X` when you run filters), or it could be as simple as the proteins having slightly different names in the two separate objects. I would recommend the first step you do is to pull out the protein data from the object in your analysis and make separate protein anndata objects, then you can filter the RNA objects and do

```python
protein_adata = protein_adata[rna_adata.obs_names]
rna_adata.obsm["protein"] = protein_adata.to_df()

```

---

<div class="post-metadata">

**Author:** ![michael101](https://avatars.discourse-cdn.com/v4/letter/m/df705f/32.png) [@michael101](https://discourse.scverse.org/u/michael101)\
**Post date:** [July 30, 2021, 3:45pm UTC](https://discourse.scverse.org/t/preparing-data-for-totalvi/149/12 "2021-07-30T15:45:59Z")

</div>

Many thanks. I think that’s actually completely right because I’ve gone back through the script and done adata.var[“feature\_types”].value\_counts() after each line and it is the sc.pp.highly\_variable\_genes step which is causing the drop.

Thanks for the advice! Much appreciated.
