{"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":"# What is about ? \n\n**Raw Data** from the competition, but denoising preprocessing by DCA package (https://scanpy.readthedocs.io/en/stable/generated/scanpy.external.pp.dca.html). See discussion: https://www.kaggle.com/competitions/open-problems-multimodal/discussion/350856\n\n**Work principle:** The main data can be accessed as f['0']['block0_values'] - without loading to RAM, but  when one takes a slice: f['0']['block0_values'][:,IX] - it is converted to numpy array (and loaded into memory?). \nThere are some details: one can do constant slices like  f['0']['block0_values'][:10,:10], but to make a  variable slice one need to process in two steps: f['0']['block0_values'][:,IX1][IX0]. \n\n**Details:** We load the data. And as an exercise and sanity check we produce plots for cell proliferation cycle visualization - so-called G1/S-G2/M plots. The results are as expected - the \"cyclic\" structes is better seen comparing as to original files (see version 3 of the present notebook). Remark: one can compare with another (more simple denoising (\"pooling\")) here: https://www.kaggle.com/code/alexandervc/mmscel-cell-cycle-03-daybydaychange-megakaryocyte#G1/S-G2/M-plots-to-analyse-the-cell-cycle\n\n**Context:** See paper: https://arxiv.org/abs/2208.05229 for cell cycle analysis. \"Computational challenges of cell cycle analysis using single cell transcriptomics\" Alexander Chervov, Andrei Zinovyev. All the work for that paper has been done on Kaggle - see datasets and code therein: https://www.kaggle.com/alexandervc/datasets?scroll=true https://www.kaggle.com/andreizinovyev Hundreds of single cell RNA seq datasets were analyzed.\n","metadata":{}},{"cell_type":"markdown","source":"# Install/import","metadata":{}},{"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","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2022-10-23T21:44:29.122044Z","iopub.execute_input":"2022-10-23T21:44:29.122430Z","iopub.status.idle":"2022-10-23T21:44:29.134990Z","shell.execute_reply.started":"2022-10-23T21:44:29.122398Z","shell.execute_reply":"2022-10-23T21:44:29.133690Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\nimport time\nt0start = time.time()\n\nimport pandas as pd\nimport numpy as np\nimport os\nimport sys\n\nimport matplotlib.pyplot as plt\n#plt.style.use('dark_background')\nimport seaborn as sns\n\n#If you see a urllib warning running this cell, go to \"Settings\" on the right hand side, \n#and turn on internet. Note, you need to be phone verified.\n!pip install --quiet tables\n\n\nimport h5py\n!pip install hdf5plugin~=2.0 # https://forum.hdfgroup.org/t/cant-open-directory-usr-local-hdf5-lib-plugin/9738/4\nimport hdf5plugin\n\n# !pip install scanpy\n# import scanpy as sc\n# import anndata\n\nDATA_DIR = \"/kaggle/input/open-problems-multimodal/\"\nFP_CELL_METADATA = os.path.join(DATA_DIR,\"metadata.csv\")\n\nFP_CITE_TRAIN_INPUTS = os.path.join(DATA_DIR,\"train_cite_inputs.h5\")\nFP_CITE_TRAIN_TARGETS = os.path.join(DATA_DIR,\"train_cite_targets.h5\")\nFP_CITE_TEST_INPUTS = os.path.join(DATA_DIR,\"test_cite_inputs.h5\")\n\nFP_MULTIOME_TRAIN_INPUTS = os.path.join(DATA_DIR,\"train_multi_inputs.h5\")\nFP_MULTIOME_TRAIN_TARGETS = os.path.join(DATA_DIR,\"train_multi_targets.h5\")\nFP_MULTIOME_TEST_INPUTS = os.path.join(DATA_DIR,\"test_multi_inputs.h5\")\n\nFP_SUBMISSION = os.path.join(DATA_DIR,\"sample_submission.csv\")\nFP_EVALUATION_IDS = os.path.join(DATA_DIR,\"evaluation_ids.csv\")\n\ndf_cell = pd.read_csv(FP_CELL_METADATA)\ndf_cell","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:45:06.557728Z","iopub.execute_input":"2022-10-23T21:45:06.558170Z","iopub.status.idle":"2022-10-23T21:45:28.421172Z","shell.execute_reply.started":"2022-10-23T21:45:06.558128Z","shell.execute_reply":"2022-10-23T21:45:28.419895Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# \"Load\" data and work with it ","metadata":{}},{"cell_type":"code","source":"print('Look on Features: ')\nfilename = '/kaggle/input/citeseq-denoised-dca/test_cite_inputs_denoised_dca.h5'\nfilename = '/kaggle/input/citeseq-denoised-dca/train_cite_inputs_denoised_dca.h5'\n#filename = '/kaggle/input/open-problems-multimodal/train_multi_inputs.h5'\n#filename = '/kaggle/input/open-problems-multimodal/train_multi_inputs.h5'\n#filename = '/kaggle/input/open-problems-multimodal/train_cite_inputs.h5'\n\n#filename = '/kaggle/input/citeseq-denoised/test_cite_inputs_denoised.h5'\n\n#filename = '/kaggle/input/open-problems-multimodal/test_cite_inputs.h5'\n\nf2 = h5py.File(filename,'r')#, mode)","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:45:43.235764Z","iopub.execute_input":"2022-10-23T21:45:43.236818Z","iopub.status.idle":"2022-10-23T21:45:43.251141Z","shell.execute_reply.started":"2022-10-23T21:45:43.236746Z","shell.execute_reply":"2022-10-23T21:45:43.249938Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(f2.keys() )","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.160524Z","iopub.status.idle":"2022-10-23T21:44:40.161152Z","shell.execute_reply.started":"2022-10-23T21:44:40.160936Z","shell.execute_reply":"2022-10-23T21:44:40.160969Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"str_data_inf = ''\nif filename == '/kaggle/input/open-problems-multimodal/train_cite_inputs.h5':\n    first_key = 'train_cite_inputs'\n    str_data_inf = ' CITEseq Train '\n# elif filename == '/kaggle/input/open-problems-multimodal/train_multi_inputs.h5':\n#     first_key = 'train_multi_inputs'\n#     str_data_inf = ' Multi Train Features '\nelif filename == '/kaggle/input/open-problems-multimodal/test_cite_inputs.h5':\n    first_key = 'test_cite_inputs'#'train_multi_inputs' \n    str_data_inf = ' CITEseq Test '\nelif filename == '/kaggle/input/citeseq-denoised-dca/test_cite_inputs_denoised_dca.h5':\n    first_key = 'test_cite_inputs_denoised_dca'\n    str_data_inf = ' MAGIC CITEseq Test '\nelif filename == '/kaggle/input/citeseq-denoised-dca/train_cite_inputs_denoised_dca.h5':\n    first_key = 'train_cite_inputs_denoised_dca'\n    str_data_inf = ' MAGIC CITEseq Train '\n# else:\n#     first_key = '0'#'train_multi_inputs' \n    \nprint( f2[first_key].keys() )","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:45:45.279394Z","iopub.execute_input":"2022-10-23T21:45:45.280214Z","iopub.status.idle":"2022-10-23T21:45:45.294655Z","shell.execute_reply.started":"2022-10-23T21:45:45.280162Z","shell.execute_reply":"2022-10-23T21:45:45.293668Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for k in f2[first_key].keys():\n    print(type(f2[first_key][k]), f2[first_key][k].shape)","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:45:49.653332Z","iopub.execute_input":"2022-10-23T21:45:49.653755Z","iopub.status.idle":"2022-10-23T21:45:49.663482Z","shell.execute_reply.started":"2022-10-23T21:45:49.653719Z","shell.execute_reply":"2022-10-23T21:45:49.662040Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"f2[first_key]['axis0'][0:10], f2[first_key]['axis1'][0:10], f2[first_key]['block0_items'][0:10],  f2[first_key]['block0_values'][0:10,0:10],\n#what's inside","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:45:50.870279Z","iopub.execute_input":"2022-10-23T21:45:50.871417Z","iopub.status.idle":"2022-10-23T21:45:50.909620Z","shell.execute_reply.started":"2022-10-23T21:45:50.871374Z","shell.execute_reply":"2022-10-23T21:45:50.908869Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"type(f2[first_key]['block0_values'][:10,:10])","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:45:53.065245Z","iopub.execute_input":"2022-10-23T21:45:53.066025Z","iopub.status.idle":"2022-10-23T21:45:53.074585Z","shell.execute_reply.started":"2022-10-23T21:45:53.065983Z","shell.execute_reply":"2022-10-23T21:45:53.073355Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print( dir(f2[first_key]['block0_values'][:10,:10]) )","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.169982Z","iopub.status.idle":"2022-10-23T21:44:40.170360Z","shell.execute_reply.started":"2022-10-23T21:44:40.170178Z","shell.execute_reply":"2022-10-23T21:44:40.170196Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"list_genes_ensembl_ids = [t.decode().split('_')[0] for t in f2[first_key]['axis0'] ]\nprint( len(list_genes_ensembl_ids), list_genes_ensembl_ids[:10] )\nlist_genes_symbols = [t.decode().split('_')[1] for t in f2[first_key]['axis0'] ]\nprint( len(list_genes_symbols), list_genes_symbols[:10] )","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:45:56.805635Z","iopub.execute_input":"2022-10-23T21:45:56.806782Z","iopub.status.idle":"2022-10-23T21:46:00.798699Z","shell.execute_reply.started":"2022-10-23T21:45:56.806732Z","shell.execute_reply":"2022-10-23T21:46:00.797468Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"list_cells_barcodes = [t.decode() for t in f2[first_key]['axis1'] ]\nprint(len(list_cells_barcodes) , list_cells_barcodes[:10] )","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:46:00.800416Z","iopub.execute_input":"2022-10-23T21:46:00.800758Z","iopub.status.idle":"2022-10-23T21:46:07.047856Z","shell.execute_reply.started":"2022-10-23T21:46:00.800725Z","shell.execute_reply":"2022-10-23T21:46:07.046636Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print( type(f2[first_key]['block0_values']), dir(f2[first_key]['block0_values'] ) ) # ???","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.178141Z","iopub.status.idle":"2022-10-23T21:44:40.178525Z","shell.execute_reply.started":"2022-10-23T21:44:40.178333Z","shell.execute_reply":"2022-10-23T21:44:40.178351Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# G1/S-G2/M cell proliferation cycle plot for all cells \n\nIt will not be very clear - we need to separate days and cell types - see next subsections.","metadata":{}},{"cell_type":"code","source":"# useful genes\n#mice elsemble\nG1S_genes_Tirosh = ['ENSMUSG00000005410', 'ENSMUSG00000027342', 'ENSMUSG00000025747', 'ENSMUSG00000024742', 'ENSMUSG00000002870', 'ENSMUSG00000022673', 'ENSMUSG00000030978', 'ENSMUSG00000029591', 'ENSMUSG00000031821', 'ENSMUSG00000026355', 'ENSMUSG00000055612', 'ENSMUSG00000037474', 'ENSMUSG00000025395', 'ENSMUSG00000001228', 'ENSMUSG00000031629', 'ENSMUSG00000025001', 'ENSMUSG00000023104', 'ENSMUSG00000028884', 'ENSMUSG00000028693', 'ENSMUSG00000030346', 'ENSMUSG00000006715', 'ENSMUSG00000027242', 'ENSMUSG00000004642', 'ENSMUSG00000028212', 'ENSMUSG00000041712', 'ENSMUSG00000030726', 'ENSMUSG00000024151', 'ENSMUSG00000022360', 'ENSMUSG00000027323', 'ENSMUSG00000020649', 'ENSMUSG00000000028', 'ENSMUSG00000017499', 'ENSMUSG00000039748', 'ENSMUSG00000032397', 'ENSMUSG00000022422', 'ENSMUSG00000030528', 'ENSMUSG00000028282', 'ENSMUSG00000028560', 'ENSMUSG00000042489', 'ENSMUSG00000006678', 'ENSMUSG00000022945', 'ENSMUSG00000034329', 'ENSMUSG00000046179']\nG2M_genes_Tirosh = ['ENSMUSG00000054717', 'ENSMUSG00000019942', 'ENSMUSG00000027306', 'ENSMUSG00000001403', 'ENSMUSG00000017716', 'ENSMUSG00000027469', 'ENSMUSG00000020914', 'ENSMUSG00000024056', 'ENSMUSG00000062248', 'ENSMUSG00000026683', 'ENSMUSG00000028044', 'ENSMUSG00000031004', 'ENSMUSG00000019961', 'ENSMUSG00000026605', 'ENSMUSG00000037313', 'ENSMUSG00000020808', 'ENSMUSG00000034349', 'ENSMUSG00000032218', 'ENSMUSG00000048327', 'ENSMUSG00000037725', 'ENSMUSG00000020897', 'ENSMUSG00000027379', 'ENSMUSG00000012443', 'ENSMUSG00000015749', 'ENSMUSG00000036752', 'ENSMUSG00000022385', 'ENSMUSG00000024795', 'ENSMUSG00000044783', 'ENSMUSG00000023505', 'ENSMUSG00000020737', 'ENSMUSG00000006398', 'ENSMUSG00000038379', 'ENSMUSG00000044201', 'ENSMUSG00000028678', 'ENSMUSG00000022391', 'ENSMUSG00000038252', 'ENSMUSG00000037544', 'ENSMUSG00000048922', 'ENSMUSG00000028873', 'ENSMUSG00000027699', 'ENSMUSG00000032254', 'ENSMUSG00000020330', 'ENSMUSG00000027496', 'ENSMUSG00000068744', 'ENSMUSG00000036777', 'ENSMUSG00000004880', 'ENSMUSG00000040549', 'ENSMUSG00000045328', 'ENSMUSG00000005698', 'ENSMUSG00000026622', 'ENSMUSG00000035293', 'ENSMUSG00000074802', 'ENSMUSG00000009575', 'ENSMUSG00000029177']\n#gene symbol\nG1S_genes_Tirosh = ['MCM5', 'PCNA', 'TYMS', 'FEN1', 'MCM2', 'MCM4', 'RRM1', 'UNG', 'GINS2', 'MCM6', 'CDCA7', 'DTL', 'PRIM1', 'UHRF1', 'MLF1IP', 'HELLS', 'RFC2', 'RPA2', 'NASP', 'RAD51AP1', 'GMNN', 'WDR76', 'SLBP', 'CCNE2', 'UBR7', 'POLD3', 'MSH2', 'ATAD2', 'RAD51', 'RRM2', 'CDC45', 'CDC6', 'EXO1', 'TIPIN', 'DSCC1', 'BLM', 'CASP8AP2', 'USP1', 'CLSPN', 'POLA1', 'CHAF1B', 'BRIP1', 'E2F8']\nG2M_genes_Tirosh = ['HMGB2', 'CDK1', 'NUSAP1', 'UBE2C', 'BIRC5', 'TPX2', 'TOP2A', 'NDC80', 'CKS2', 'NUF2', 'CKS1B', 'MKI67', 'TMPO', 'CENPF', 'TACC3', 'FAM64A', 'SMC4', 'CCNB2', 'CKAP2L', 'CKAP2', 'AURKB', 'BUB1', 'KIF11', 'ANP32E', 'TUBB4B', 'GTSE1', 'KIF20B', 'HJURP', 'CDCA3', 'HN1', 'CDC20', 'TTK', 'CDC25C', 'KIF2C', 'RANGAP1', 'NCAPD2', 'DLGAP5', 'CDCA2', 'CDCA8', 'ECT2', 'KIF23', 'HMMR', 'AURKA', 'PSRC1', 'ANLN', 'LBR', 'CKAP5', 'CENPE', 'CTCF', 'NEK2', 'G2E3', 'GAS2L3', 'CBX5', 'CENPA']\nlist_genes_fastCCsign = ['CDK1', 'UBE2C', 'TOP2A', 'TMPO', 'HJURP', 'RRM1', 'RAD51AP1', 'RRM2', 'CDC45', 'BLM', 'BRIP1', 'E2F8', 'HIST2H2AC']\n#human ensemble\nG1S_genes_Tirosh = ['ENSG00000100297', 'ENSG00000132646', 'ENSG00000176890', 'ENSG00000168496', 'ENSG00000073111', 'ENSG00000104738', 'ENSG00000167325', 'ENSG00000076248', 'ENSG00000131153', 'ENSG00000076003', 'ENSG00000144354', 'ENSG00000143476', 'ENSG00000198056', 'ENSG00000276043', 'ENSG00000151725', 'ENSG00000119969', 'ENSG00000049541', 'ENSG00000117748', 'ENSG00000132780', 'ENSG00000111247', 'ENSG00000112312', 'ENSG00000092470', 'ENSG00000163950', 'ENSG00000175305', 'ENSG00000012963', 'ENSG00000077514', 'ENSG00000095002', 'ENSG00000156802', 'ENSG00000051180', 'ENSG00000171848', 'ENSG00000093009', 'ENSG00000094804', 'ENSG00000174371', 'ENSG00000075131', 'ENSG00000136982', 'ENSG00000197299', 'ENSG00000118412', 'ENSG00000162607', 'ENSG00000092853', 'ENSG00000101868', 'ENSG00000159259', 'ENSG00000136492', 'ENSG00000129173']\nG2M_genes_Tirosh = ['ENSG00000164104', 'ENSG00000170312', 'ENSG00000137804', 'ENSG00000175063', 'ENSG00000089685', 'ENSG00000088325', 'ENSG00000131747', 'ENSG00000080986', 'ENSG00000123975', 'ENSG00000143228', 'ENSG00000173207', 'ENSG00000148773', 'ENSG00000120802', 'ENSG00000117724', 'ENSG00000013810', 'ENSG00000129195', 'ENSG00000113810', 'ENSG00000157456', 'ENSG00000169607', 'ENSG00000136108', 'ENSG00000178999', 'ENSG00000169679', 'ENSG00000138160', 'ENSG00000143401', 'ENSG00000188229', 'ENSG00000075218', 'ENSG00000138182', 'ENSG00000123485', 'ENSG00000111665', 'ENSG00000189159', 'ENSG00000117399', 'ENSG00000112742', 'ENSG00000158402', 'ENSG00000142945', 'ENSG00000100401', 'ENSG00000010292', 'ENSG00000126787', 'ENSG00000184661', 'ENSG00000134690', 'ENSG00000114346', 'ENSG00000137807', 'ENSG00000072571', 'ENSG00000087586', 'ENSG00000134222', 'ENSG00000011426', 'ENSG00000143815', 'ENSG00000175216', 'ENSG00000138778', 'ENSG00000102974', 'ENSG00000117650', 'ENSG00000092140', 'ENSG00000139354', 'ENSG00000094916', 'ENSG00000115163']\nlist_genes_fastCCsign = ['ENSG00000170312', 'ENSG00000175063', 'ENSG00000131747', 'ENSG00000120802', 'ENSG00000123485', 'ENSG00000167325', 'ENSG00000111247', 'ENSG00000171848', 'ENSG00000093009', 'ENSG00000197299', 'ENSG00000136492', 'ENSG00000129173', 'ENSG00000184260']\n\n\n","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:46:07.050061Z","iopub.execute_input":"2022-10-23T21:46:07.050491Z","iopub.status.idle":"2022-10-23T21:46:07.071321Z","shell.execute_reply.started":"2022-10-23T21:46:07.050448Z","shell.execute_reply":"2022-10-23T21:46:07.070500Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#print(G1S_genes_Tirosh)\nIX = np.where(pd.Series(list_genes_ensembl_ids ).isin(G1S_genes_Tirosh) >0) [0]\nprint(IX)\nlen(IX), len( G1S_genes_Tirosh)","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:46:07.342571Z","iopub.execute_input":"2022-10-23T21:46:07.343700Z","iopub.status.idle":"2022-10-23T21:46:07.358898Z","shell.execute_reply.started":"2022-10-23T21:46:07.343647Z","shell.execute_reply":"2022-10-23T21:46:07.357697Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nv1 = f2[first_key]['block0_values'][:,IX].sum(axis = 1 )","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.183505Z","iopub.status.idle":"2022-10-23T21:44:40.183906Z","shell.execute_reply.started":"2022-10-23T21:44:40.183690Z","shell.execute_reply":"2022-10-23T21:44:40.183708Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"IX = np.where(pd.Series(list_genes_ensembl_ids ).isin(G2M_genes_Tirosh) >0) [0]\nlen(IX), len( G2M_genes_Tirosh)","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.184979Z","iopub.status.idle":"2022-10-23T21:44:40.185320Z","shell.execute_reply.started":"2022-10-23T21:44:40.185147Z","shell.execute_reply":"2022-10-23T21:44:40.185163Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nv2 = f2[first_key]['block0_values'][:,IX].sum(axis = 1 )\nprint( len(v2), type(v2), v2.shape) ","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.186689Z","iopub.status.idle":"2022-10-23T21:44:40.187083Z","shell.execute_reply.started":"2022-10-23T21:44:40.186891Z","shell.execute_reply":"2022-10-23T21:44:40.186910Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sns.scatterplot(x = v1, y=v2)\nplt.title(str_data_inf)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.188197Z","iopub.status.idle":"2022-10-23T21:44:40.188572Z","shell.execute_reply.started":"2022-10-23T21:44:40.188377Z","shell.execute_reply":"2022-10-23T21:44:40.188395Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# G1/S-G2/M plot for day=1, cell_type = MkP (Megakaryocyte Progenitors)","metadata":{}},{"cell_type":"code","source":"#need to build a table with desired genes from input\n# values - \n#genes FOXO1, CCNE1, GAPDH, DRD4\ntest_gene_set = ['CCNE1', 'FOXO1', 'GAPDH', 'DRD4']\nIX = np.where(pd.Series(list_genes_symbols ).isin(test_gene_set) >0) [0]\ndf_cite_train_y = pd.DataFrame(f2[first_key]['block0_values'][:,IX], columns = test_gene_set)\ndf_cells_barcodes = pd.Series(list_cells_barcodes, name = 'cell_id')\ndf_cite_train_y = pd.concat([df_cells_barcodes.reset_index(drop=True), df_cite_train_y.reset_index(drop=True)], axis=1)\n\n\n#genes = list(df_cite_train_y.keys())\ndf_cell_hue = df_cell.merge(df_cite_train_y, on = 'cell_id')\ndisplay(df_cell_hue.head(15))","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:46:13.029184Z","iopub.execute_input":"2022-10-23T21:46:13.029599Z","iopub.status.idle":"2022-10-23T21:47:34.580555Z","shell.execute_reply.started":"2022-10-23T21:46:13.029566Z","shell.execute_reply":"2022-10-23T21:47:34.578687Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"df_cite_train_y = pd.read_hdf('../input/open-problems-multimodal/train_cite_targets.h5')\ngenes = list(df_cite_train_y.keys())\ndf_cell_hue = df_cell.merge(df_cite_train_y, on = 'cell_id')\ndisplay(df_cell_hue.head(15))\n#print(genes)","metadata":{"execution":{"iopub.status.busy":"2022-10-23T20:35:42.710283Z","iopub.execute_input":"2022-10-23T20:35:42.710980Z","iopub.status.idle":"2022-10-23T20:35:43.682059Z","shell.execute_reply.started":"2022-10-23T20:35:42.710935Z","shell.execute_reply":"2022-10-23T20:35:43.680945Z"}}},{"cell_type":"code","source":"type(df_cite_train_y.keys())","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:47:34.583693Z","iopub.execute_input":"2022-10-23T21:47:34.584065Z","iopub.status.idle":"2022-10-23T21:47:34.593287Z","shell.execute_reply.started":"2022-10-23T21:47:34.584031Z","shell.execute_reply":"2022-10-23T21:47:34.591968Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# G1/S-G2/M plot for day=2, cell_type = MkP (Megakaryocyte Progenitors)\n\nCompare with: https://www.kaggle.com/code/alexandervc/mmscel-cell-cycle-03-daybydaychange-megakaryocyte?scriptVersionId=105307517&cellId=12\n\n","metadata":{}},{"cell_type":"code","source":"palette = 'rainbow'#'viridis'\nn_x_subplots = 4\nday = 2\n\n# filtering only day 2 & cell type\nmask = (df_cell_hue['day'] == day) & (df_cell_hue['cell_type'] == 'MkP' )                      \ncell_ids = list(df_cell_hue['cell_id'][mask])\n\n#generating indexes for genes from defined list in G1/S\nIX = np.where(pd.Series(list_genes_ensembl_ids).isin(G1S_genes_Tirosh) >0) [0]\nv1 = (f2[first_key]['block0_values'][:,IX]).sum(axis = 1 )\n# Remark: Error - if write [IX0,IX] - \"TypeError: Only one indexing vector or array is currently allowed for fancy indexing\"\n# so we proceed as above - in two steps\n\n# same for genes in G2/M\nIX = np.where(pd.Series(list_genes_ensembl_ids ).isin(G2M_genes_Tirosh) >0) [0]\nv2 = (f2[first_key]['block0_values'][:,IX]).sum(axis = 1 )\n\nIX0 = np.where( pd.Series(list_cells_barcodes).isin(cell_ids) > 0)[0]\n\n\nc = 0\nfor gene_id in test_gene_set: # ['day', 'donor', 'cell_type']:# , 'technology']:\n    if c % n_x_subplots == 0:\n        if c > 0:\n            plt.show()\n        fig = plt.figure(figsize = (20,5) ); c = 0\n        plt.xlabel('G1/S score', fontsize = 20 )\n        plt.ylabel('G2/M score', fontsize = 20 )\n\n    c += 1\n    fig.add_subplot(1,n_x_subplots ,c)    \n\n    v = df_cell_hue[gene_id]\n\n    ax = sns.scatterplot(x=v1[IX0], y=v2[IX0], hue=v[IX0], palette = palette )\n    #plt.setp(ax.get_legend().get_texts(), fontsize=10) # for legend text\n    #plt.setp(ax.get_legend().get_title(), fontsize=10) # for legend title    \n    plt.title('Colored by '+gene_id, fontsize = 10)\n\n#plt.title(str_data_inf + 'Day 2, MkP (Megakaryocyte Progenitors) ', fontsize = 20)\n\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:47:34.594660Z","iopub.execute_input":"2022-10-23T21:47:34.595889Z","iopub.status.idle":"2022-10-23T21:48:03.510778Z","shell.execute_reply.started":"2022-10-23T21:47:34.595851Z","shell.execute_reply":"2022-10-23T21:48:03.509880Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# G1/S-G2/M plot for day=3, cell_type = MkP (Megakaryocyte Progenitors)","metadata":{}},{"cell_type":"code","source":"palette = 'rainbow'#'viridis'\nn_x_subplots = 4\nday = 3\n\nmask = (df_cell_hue['day'] == day) & (df_cell_hue['cell_type'] == 'MkP' )                      \ncell_ids = list(df_cell_hue['cell_id'][mask])\n\nIX = np.where(pd.Series(list_genes_ensembl_ids ).isin(G1S_genes_Tirosh) >0) [0]\nv1 = (f2[first_key]['block0_values'][:,IX]).sum(axis = 1 )\n# Remark: Error - if write [IX0,IX] - \"TypeError: Only one indexing vector or array is currently allowed for fancy indexing\"\n# so we proceed as above - in two steps\n\nIX = np.where(pd.Series(list_genes_ensembl_ids ).isin(G2M_genes_Tirosh) >0) [0]\nv2 = (f2[first_key]['block0_values'][:,IX]).sum(axis = 1 )\n\nIX0 = np.where( pd.Series(list_cells_barcodes).isin(cell_ids) > 0)[0]\n\n\nc = 0\nfor gene_id in  test_gene_set: #genes: # ['day', 'donor', 'cell_type']:# , 'technology']:\n    if c % n_x_subplots == 0:\n        if c > 0:\n            plt.show()\n        fig = plt.figure(figsize = (20,5) ); c = 0\n        plt.xlabel('G1/S score', fontsize = 20 )\n        plt.ylabel('G2/M score', fontsize = 20 )\n\n    c += 1\n    fig.add_subplot(1,n_x_subplots ,c)    \n\n    v = df_cell_hue[gene_id]\n\n    ax = sns.scatterplot(x=v1[IX0], y=v2[IX0], hue=v[IX0], palette = palette )\n    plt.setp(ax.get_legend().get_texts(), fontsize=10) # for legend text\n    plt.setp(ax.get_legend().get_title(), fontsize=10) # for legend title    \n    plt.title('Colored by '+gene_id, fontsize = 10)\n\n#plt.title(str_data_inf + 'Day 4, MkP (Megakaryocyte Progenitors) ', fontsize = 20)\n\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.195231Z","iopub.status.idle":"2022-10-23T21:44:40.195622Z","shell.execute_reply.started":"2022-10-23T21:44:40.195431Z","shell.execute_reply":"2022-10-23T21:44:40.195449Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# G1/S-G2/M plot for day=4, cell_type = MkP (Megakaryocyte Progenitors)","metadata":{}},{"cell_type":"code","source":"palette = 'rainbow'#'viridis'\nn_x_subplots = 4\nday = 4\n\nmask = (df_cell_hue['day'] == day) & (df_cell_hue['cell_type'] == 'MkP' )                      \ncell_ids = list(df_cell_hue['cell_id'][mask])\n\nIX = np.where(pd.Series(list_genes_ensembl_ids ).isin(G1S_genes_Tirosh) >0) [0]\nv1 = (f2[first_key]['block0_values'][:,IX]).sum(axis = 1 )\n# Remark: Error - if write [IX0,IX] - \"TypeError: Only one indexing vector or array is currently allowed for fancy indexing\"\n# so we proceed as above - in two steps\n\nIX = np.where(pd.Series(list_genes_ensembl_ids ).isin(G2M_genes_Tirosh) >0) [0]\nv2 = (f2[first_key]['block0_values'][:,IX]).sum(axis = 1 )\n\nIX0 = np.where( pd.Series(list_cells_barcodes).isin(cell_ids) > 0)[0]\n\n\nc = 0\nfor gene_id in test_gene_set: #genes: # ['day', 'donor', 'cell_type']:# , 'technology']:\n    if c % n_x_subplots == 0:\n        if c > 0:\n            plt.show()\n        fig = plt.figure(figsize = (20,5) ); c = 0\n        plt.xlabel('G1/S score', fontsize = 20 )\n        plt.ylabel('G2/M score', fontsize = 20 )\n\n    c += 1\n    fig.add_subplot(1,n_x_subplots ,c)    \n\n    v = df_cell_hue[gene_id]\n\n    ax = sns.scatterplot(x=v1[IX0], y=v2[IX0], hue=v[IX0], palette = palette )\n    plt.setp(ax.get_legend().get_texts(), fontsize=10) # for legend text\n    plt.setp(ax.get_legend().get_title(), fontsize=10) # for legend title    \n    plt.title('Colored by '+gene_id, fontsize = 10)\n\n#plt.title(str_data_inf + 'Day 4, MkP (Megakaryocyte Progenitors) ', fontsize = 20)\n\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.196721Z","iopub.status.idle":"2022-10-23T21:44:40.197134Z","shell.execute_reply.started":"2022-10-23T21:44:40.196940Z","shell.execute_reply":"2022-10-23T21:44:40.196965Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Loop over all cell types ","metadata":{}},{"cell_type":"code","source":"dict_cell_types = { 'MasP': 'Mast Cell Progenitor',\n'MkP': 'Megakaryocyte Progenitor',\n'NeuP': 'Neutrophil Progenitor',\n'MoP': 'Monocyte Progenitor',\n'EryP': 'Erythrocyte Progenitor',\n'HSC': 'Hematoploetic Stem Cell',\n'BP': 'B-Cell Progenitor'}","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.199169Z","iopub.status.idle":"2022-10-23T21:44:40.199560Z","shell.execute_reply.started":"2022-10-23T21:44:40.199363Z","shell.execute_reply":"2022-10-23T21:44:40.199382Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_cell.head(1)\nd = pd.DataFrame(index =list_cells_barcodes,  )\nd = d.join(df_cell.set_index('cell_id'), how = 'left')\nd['str_donor'] = d['donor'].apply(lambda x: str(x))\nd","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.200725Z","iopub.status.idle":"2022-10-23T21:44:40.201126Z","shell.execute_reply.started":"2022-10-23T21:44:40.200941Z","shell.execute_reply":"2022-10-23T21:44:40.200959Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for cell_type in ['MkP','HSC',  'EryP', 'NeuP', 'MasP' ,  'MoP' , 'BP' ]:\n    print();print();print()\n    print(cell_type)\n    \n    for day in [2,3,4,7]:\n\n        mask = (df_cell['day'] == day) & (df_cell['cell_type'] == cell_type )                      \n        #print('Day',day, mask.sum() )\n        cell_ids = list(df_cell['cell_id'][mask])\n        #print(len(cell_ids) , cell_ids[:10])\n\n        IX0 = np.where( pd.Series(list_cells_barcodes).isin(cell_ids) > 0)[0]\n        print('Day', day,mask.sum(), len(IX0))# , IX0[:10]) # , len(cell_ids ))\n\n        plt.figure(figsize = (20,10))\n        if d['str_donor'].iloc[IX0].nunique() > 1:\n            sns.scatterplot(x = v1[IX0], y=v2[IX0] , hue = d['str_donor'].iloc[IX0], palette = 'rainbow')#  ['blue'])\n        else:\n            sns.scatterplot(x = v1[IX0], y=v2[IX0] , hue = d['str_donor'].iloc[IX0], palette =  ['blue'])\n            \n        plt.xlabel('G1/S score', fontsize = 20 )\n        plt.ylabel('G2/M score', fontsize = 20 )\n        plt.title(str_data_inf + '  Day '+str(day) + ', '+ cell_type +  ' ' + dict_cell_types[cell_type] + ' n_cells=%d'%len(IX0) , fontsize = 20)\n\n        plt.show()","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.202425Z","iopub.status.idle":"2022-10-23T21:44:40.202816Z","shell.execute_reply.started":"2022-10-23T21:44:40.202619Z","shell.execute_reply":"2022-10-23T21:44:40.202637Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('%.1f seconds passed total '%(time.time()-t0start) )\n","metadata":{"execution":{"iopub.status.busy":"2022-10-23T21:44:40.204268Z","iopub.status.idle":"2022-10-23T21:44:40.204643Z","shell.execute_reply.started":"2022-10-23T21:44:40.204457Z","shell.execute_reply":"2022-10-23T21:44:40.204475Z"},"trusted":true},"execution_count":null,"outputs":[]}]}