{"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":"# TOC \n### (This notebook is under construction !!!)\n\n---\n\n#### Version 1\n##### - parquet -> adata\n##### - pre-preprocessing\n##### - UMAP\n\n#### Version 2\n##### - Add filter based on `adata.var.highly_variable_nbatches`\n\n#### Version 3\n##### - Add `sc.pl.embedding_density`","metadata":{"execution":{"iopub.status.busy":"2023-09-20T07:12:14.963152Z","iopub.execute_input":"2023-09-20T07:12:14.963702Z","iopub.status.idle":"2023-09-20T07:12:14.970688Z","shell.execute_reply.started":"2023-09-20T07:12:14.963663Z","shell.execute_reply":"2023-09-20T07:12:14.969423Z"}}},{"cell_type":"code","source":"!pip install scanpy","metadata":{"execution":{"iopub.status.busy":"2023-09-22T16:46:06.533477Z","iopub.execute_input":"2023-09-22T16:46:06.533917Z","iopub.status.idle":"2023-09-22T16:46:27.590796Z","shell.execute_reply.started":"2023-09-22T16:46:06.533887Z","shell.execute_reply":"2023-09-22T16:46:27.588822Z"},"_kg_hide-output":true,"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python Docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load\n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\n\n# Input data files are available in the read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\nimport os\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    for filename in filenames:\n        print(os.path.join(dirname, filename))\n\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session\n\nimport gc\nimport scanpy as sc\n\nimport matplotlib.pyplot as plt","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2023-09-22T16:46:27.594703Z","iopub.execute_input":"2023-09-22T16:46:27.595103Z","iopub.status.idle":"2023-09-22T16:46:35.305540Z","shell.execute_reply.started":"2023-09-22T16:46:27.595066Z","shell.execute_reply":"2023-09-22T16:46:35.303849Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"---","metadata":{}},{"cell_type":"markdown","source":"# parquet -> adata\n### run on local machine (~1h), not have enough memory on kaggle kernel\n","metadata":{}},{"cell_type":"markdown","source":"```python\nimport scanpy as sc\n\nadata_obs_meta = pd.read_csv('/kaggle/input/open-problems-single-cell-perturbations/adata_obs_meta.csv').set_index('obs_id')\nadata_train = pd.read_parquet('/kaggle/input/open-problems-single-cell-perturbations/adata_train.parquet')\n\n# parquet to matrix\nadata_xdf = adata_train.pivot_table(index=\"obs_id\", columns=\"gene\", values=\"count\").fillna(0)\n\n# matrix to adata\nadata = sc.AnnData(\n    X=adata_xdf,\n    obs=adata_obs_meta.loc[adata_xdf.index],\n    var=pd.DataFrame(adata_xdf.columns, index=adata_xdf.columns)\n)\n\n# Save h5ad (~40G)\nadata.write_h5ad('adata.h5ad')\n```","metadata":{"execution":{"iopub.status.busy":"2023-09-19T06:56:55.997787Z","iopub.execute_input":"2023-09-19T06:56:55.999354Z","iopub.status.idle":"2023-09-19T06:57:18.035646Z","shell.execute_reply.started":"2023-09-19T06:56:55.999306Z","shell.execute_reply":"2023-09-19T06:57:18.034329Z"}}},{"cell_type":"markdown","source":"---","metadata":{}},{"cell_type":"markdown","source":"# pre-preprocessing\n### Please help if you are an expert on scRNA-seq !!!\n### I'm not sure if this is a good (standerd) pre-preprocessing pipeline.","metadata":{}},{"cell_type":"markdown","source":"```python\n\n# Normalize each cell by total counts over all genes, inplace\n# ref: https://github.com/openproblems-bio/openproblems-v2/\nsc.pp.normalize_total(adata=adata, target_sum=1e6, inplace=True)\nsc.pp.log1p(adata)\n\n# HVG\nsc.pp.highly_variable_genes(adata, n_top_genes=1000, batch_key='donor_id')\n\n# PCA\nsc.tl.pca(adata, n_comps=50, use_highly_variable=True)\n\n# UMAP\nsc.pp.neighbors(adata, n_neighbors=10, n_pcs=20)\nsc.tl.umap(adata)\n\n# Select HVG\nhvg_adata = adata[:, adata.var.highly_variable]\n\n# Save h5ad (~2G)\nhvg_adata.write_h5ad('hvg.h5ad')\n```","metadata":{"execution":{"iopub.status.busy":"2023-09-20T06:42:27.465476Z","iopub.execute_input":"2023-09-20T06:42:27.465901Z","iopub.status.idle":"2023-09-20T06:42:27.480183Z","shell.execute_reply.started":"2023-09-20T06:42:27.465871Z","shell.execute_reply":"2023-09-20T06:42:27.478292Z"}}},{"cell_type":"markdown","source":"---","metadata":{}},{"cell_type":"code","source":"hvg_adata = sc.read('/kaggle/input/op2-adata-train-h5ad/hvg.h5ad')\nhvg_adata.obs['control_tag'] = hvg_adata.obs['control'].astype(str)\nhvg_adata.obs['col_tag'] = hvg_adata.obs['col'].astype(str)","metadata":{"execution":{"iopub.status.busy":"2023-09-22T16:46:35.306868Z","iopub.execute_input":"2023-09-22T16:46:35.307689Z","iopub.status.idle":"2023-09-22T16:47:04.436946Z","shell.execute_reply.started":"2023-09-22T16:46:35.307645Z","shell.execute_reply":"2023-09-22T16:47:04.432329Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"hvg_adata","metadata":{"execution":{"iopub.status.busy":"2023-09-22T16:47:04.452380Z","iopub.execute_input":"2023-09-22T16:47:04.454068Z","iopub.status.idle":"2023-09-22T16:47:04.496006Z","shell.execute_reply.started":"2023-09-22T16:47:04.453864Z","shell.execute_reply":"2023-09-22T16:47:04.493897Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# ref: https://www.kaggle.com/competitions/open-problems-single-cell-perturbations/discussion/440956#2449975\n# hvg_adata = hvg_adata[:, hvg_adata.var.highly_variable_nbatches>=2].copy()\nhvg_adata = hvg_adata[:, hvg_adata.var.highly_variable_nbatches>=3].copy()","metadata":{"execution":{"iopub.status.busy":"2023-09-22T16:47:04.499609Z","iopub.execute_input":"2023-09-22T16:47:04.500148Z","iopub.status.idle":"2023-09-22T16:47:15.260608Z","shell.execute_reply.started":"2023-09-22T16:47:04.500110Z","shell.execute_reply":"2023-09-22T16:47:15.258614Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# PCA\nsc.tl.pca(hvg_adata, n_comps=50, use_highly_variable=True)\n\n# UMAP\nsc.pp.neighbors(hvg_adata, n_neighbors=10, n_pcs=20)\nsc.tl.umap(hvg_adata)","metadata":{"execution":{"iopub.status.busy":"2023-09-22T16:47:15.262967Z","iopub.execute_input":"2023-09-22T16:47:15.263371Z","iopub.status.idle":"2023-09-22T16:54:50.353766Z","shell.execute_reply.started":"2023-09-22T16:47:15.263338Z","shell.execute_reply":"2023-09-22T16:54:50.351537Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sc.pl.umap(hvg_adata, color='cell_type')","metadata":{"execution":{"iopub.status.busy":"2023-09-22T16:54:50.356025Z","iopub.execute_input":"2023-09-22T16:54:50.357283Z","iopub.status.idle":"2023-09-22T16:54:52.652060Z","shell.execute_reply.started":"2023-09-22T16:54:50.357241Z","shell.execute_reply":"2023-09-22T16:54:52.650487Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sc.pl.umap(hvg_adata, color='donor_id')","metadata":{"execution":{"iopub.status.busy":"2023-09-22T16:54:52.654480Z","iopub.execute_input":"2023-09-22T16:54:52.654979Z","iopub.status.idle":"2023-09-22T16:54:54.707370Z","shell.execute_reply.started":"2023-09-22T16:54:52.654937Z","shell.execute_reply":"2023-09-22T16:54:54.705669Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sc.pl.umap(hvg_adata, color='plate_name')","metadata":{"execution":{"iopub.status.busy":"2023-09-22T16:54:54.709478Z","iopub.execute_input":"2023-09-22T16:54:54.710273Z","iopub.status.idle":"2023-09-22T16:54:56.883489Z","shell.execute_reply.started":"2023-09-22T16:54:54.710225Z","shell.execute_reply":"2023-09-22T16:54:56.882057Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sc.pl.umap(hvg_adata, color='control_tag')","metadata":{"execution":{"iopub.status.busy":"2023-09-22T16:54:56.889649Z","iopub.execute_input":"2023-09-22T16:54:56.890205Z","iopub.status.idle":"2023-09-22T16:54:58.857984Z","shell.execute_reply.started":"2023-09-22T16:54:56.890160Z","shell.execute_reply":"2023-09-22T16:54:58.856651Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sc.pl.umap(hvg_adata, color='row')","metadata":{"execution":{"iopub.status.busy":"2023-09-22T16:54:58.860211Z","iopub.execute_input":"2023-09-22T16:54:58.860670Z","iopub.status.idle":"2023-09-22T16:55:01.021059Z","shell.execute_reply.started":"2023-09-22T16:54:58.860623Z","shell.execute_reply":"2023-09-22T16:55:01.019400Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sc.pl.umap(hvg_adata, color='col_tag')","metadata":{"execution":{"iopub.status.busy":"2023-09-22T16:55:01.023126Z","iopub.execute_input":"2023-09-22T16:55:01.023553Z","iopub.status.idle":"2023-09-22T16:55:03.339431Z","shell.execute_reply.started":"2023-09-22T16:55:01.023520Z","shell.execute_reply":"2023-09-22T16:55:03.338153Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"---","metadata":{}},{"cell_type":"code","source":"sc.tl.embedding_density(hvg_adata, groupby='cell_type')\nsc.tl.embedding_density(hvg_adata, groupby='donor_id')\nsc.tl.embedding_density(hvg_adata, groupby='plate_name')\nsc.tl.embedding_density(hvg_adata, groupby='control_tag')","metadata":{"execution":{"iopub.status.busy":"2023-09-22T16:55:03.341147Z","iopub.execute_input":"2023-09-22T16:55:03.341539Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sc.pl.embedding_density(hvg_adata, groupby='cell_type', ncols=3)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sc.pl.embedding_density(hvg_adata, groupby='donor_id', ncols=3)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sc.pl.embedding_density(hvg_adata, groupby='plate_name', ncols=2)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sc.pl.embedding_density(hvg_adata, groupby='control_tag', ncols=2)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"---","metadata":{}},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}