{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Integration of multi-modal data using MultiVI\n\nMultiVI `[Ashuach et al., 2021]` is used to integrate multiomic datasets with single-modality (expression or accessibility) datasets. We are going to generate a latent embedding of the multiomic datasets for each cell type.","metadata":{}},{"cell_type":"markdown","source":"## Setting up environment","metadata":{}},{"cell_type":"code","source":"!pip install scvi-tools scanpy","metadata":{"execution":{"iopub.status.busy":"2023-10-25T06:35:21.426882Z","iopub.execute_input":"2023-10-25T06:35:21.427936Z","iopub.status.idle":"2023-10-25T06:36:11.544552Z","shell.execute_reply.started":"2023-10-25T06:35:21.427901Z","shell.execute_reply":"2023-10-25T06:36:11.543118Z"},"collapsed":true,"jupyter":{"outputs_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import os\nos.environ['CUDA_VISIBLE_DEVICES'] = '0'","metadata":{"execution":{"iopub.status.busy":"2023-10-25T06:37:54.609302Z","iopub.execute_input":"2023-10-25T06:37:54.610439Z","iopub.status.idle":"2023-10-25T06:37:54.616277Z","shell.execute_reply.started":"2023-10-25T06:37:54.610400Z","shell.execute_reply":"2023-10-25T06:37:54.615044Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import scanpy as sc\nimport numpy as np\nimport pandas as pd\nimport scvi\nimport scvi.model\nimport torch\n\ntorch.set_float32_matmul_precision('medium')","metadata":{"execution":{"iopub.status.busy":"2023-10-25T06:36:29.529176Z","iopub.execute_input":"2023-10-25T06:36:29.529603Z","iopub.status.idle":"2023-10-25T06:36:56.062396Z","shell.execute_reply.started":"2023-10-25T06:36:29.529569Z","shell.execute_reply":"2023-10-25T06:36:56.061212Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Loading dataset\n\nThe preprocessed Multiome dataset with `h5ad` extension used for saving/loading anndata objects can be loaded using the `sc.read_h5ad` function. The datasets is not normalized or preprocessed.","metadata":{}},{"cell_type":"code","source":"data_dir = '/kaggle/input/open-problems-single-cell-perturbations'\nadata_mvi = sc.read_h5ad(os.path.join(data_dir, 'multiome_train.h5ad'))\nn_genes = sum(adata_mvi.var['feature_type'] == 'Gene Expression')\nn_regions = sum(adata_mvi.var['feature_type'] == 'Peaks')\nadata_mvi.X = adata_mvi.X.toarray()\n\nprint('n_genes: ', n_genes)\nprint('n_regions: ', n_regions)\nprint(adata_mvi)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Dataset setup\nNow is the time to setup the dataset for MultiVI to prepare the dataset for training. Just like scvi-tools models, you can call `scvi.model.MULTIVI.setup_anndata` to setup your data. This function will accept the following arguments:\n\n- `adata`: AnnData object containing the data to be preprocessed\n- `batch_key`: The key in `adata.obs` that contains the batch information\n- `categorical_covariate_keys`: A list of keys in `adata.obs` that contain categorical covariates\n- `continuous_covariate_keys`: A list of keys in `adata.obs` that contain continuous covariates\n- `protein_expression_obsm_key`: key in `adata.obsm` for protein expression data.\n- `protein_names_uns_key`: key in `adata.uns` for protein names. If None, will use the column names of `adata.obsm[protein_expression_obsm_key]` if it is a DataFrame, else will assign sequential names to proteins.","metadata":{}},{"cell_type":"code","source":"scvi.model.MULTIVI.setup_anndata(adata_mvi, batch_key='donor_id')","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## MultiVI Model\n\nYou can create a MultiVI model by creating an object from `scvi.model.MULTIVI` class. The constructor of this class takes the following arguments:\n\n**Data related parameters:** \n- `adata`: AnnData object containing train/valid/test data\n\n**Model architecture parameters:**\n- `n_genes`: Number of genes in the expression data\n- `n_regions`: Number of regions in the accessibility data \n- `n_latent`: Dimensionality of the latent space. If `None`, defaults to square root\n        of `n_hidden`.\n- `n_hidden`: Number of nodes per hidden layer. If `None`, defaults to square root\n        of number of regions.\n- `modality_weight`: Weighting scheme across modalities. One of the following:\n   - ``\"equal\"``: Equal weight in each modality\n   - ``\"universal\"``: Learn weights across modalities w_m.\n   - ``\"cell\"``: Learn weights across modalities and cells. w_{m,c}\n- `modality_penalty`: Training Penalty across modalities. One of the following:\n   - ``\"Jeffreys\"``: Jeffreys penalty to align modalities\n   - ``\"MMD\"``: MMD penalty to align modalities\n   - ``\"None\"``: No penalty\n- `dropout_rate`: Dropout rate for neural networks.\n- `latent_distribution`: Distribution of the latent space. One of the following:\n   - ``\"normal\"``: Normal distribution\n   - ``\"ln\"``: Log-normal distribution\n- `use_batch_norm`: Whether to use batch norm in the neural networks. Only of the following:\n   - ``\"none\"``: No batch norm\n   - ``\"encoder\"``: Batch norm in the encoder\n   - ``\"decoder\"``: Batch norm in the decoder\n   - ``\"both\"``: Batch norm in both encoder and decoder\n- `use_layer_norm`: Whether to use layer norm in the neural networks. Only of the following:\n   - ``\"none\"``: No layer norm\n   - ``\"encoder\"``: Layer norm in the encoder\n   - ``\"decoder\"``: Layer norm in the decoder\n   - ``\"both\"``: Layer norm in both encoder and decoder","metadata":{}},{"cell_type":"code","source":"model_params = dict(\n    n_genes = n_genes,\n    n_regions = n_regions,\n    n_latent = 64,\n    n_hidden = 512,\n    modality_weights = 'cell',\n    modality_penalty = 'Jeffreys',\n    gene_likelihood = 'zinb',\n    use_batch_norm = 'none',\n    use_layer_norm = 'both',\n    latent_distribution = 'normal',\n    dropout_rate = 0.1,\n)\nmvi = scvi.model.MULTIVI(adata_mvi, **model_params)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Training MultiVI\n\nIn order to train your MultiVI model, you need to use `train` function of your `model`. This function accepts the following parameters:\n- `max_epochs`: Maximum number of epochs to train the model. CPA generally converges after high number of epochs, so you can set this to a high value.\n- `use_gpu`: If you have a GPU, you can set this to `True` to speed up the training process.\n- `batch_size`: Batch size for training. You can set this to a high value (e.g. 512, 1024, 2048) if you have a GPU. \n- `lr`: Learning rate for training.\n- `weight_decay`: Weight decay for training.\n- `early_stopping`: Whether to use early stopping or not.\n- `early_stopping_patience`: Number of epochs to wait before stopping training if validation metric does not improve.\n- `check_val_every_n_epoch`: Number of epochs to wait before running validation.\n","metadata":{}},{"cell_type":"code","source":"train_params = dict(\n    max_epochs = 500,\n    lr = 1e-3,\n    use_gpu = True,\n    batch_size = 8,\n    weight_decay = 1e-5,\n    check_val_every_n_epoch = 1,\n    early_stopping = True,\n    save_best = True,\n)\nmvi.train(**train_params)\nmvi.save('models/mvi', overwrite=True)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"mvi: scvi.model.MULTIVI = scvi.model.MULTIVI.load('models/mvi', adata_mvi)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"latent_outputs = mvi.get_latent_representation()","metadata":{},"execution_count":null,"outputs":[]}]}