{"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\nEDA for multiome targets (scRNA-seq) data.\nCell cycle G1/S-G2/M plots (kind of sanity check that data fit some  expectated patterns). \n\n\n\nSee 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.\nSee also discussion https://www.kaggle.com/competitions/open-problems-multimodal/discussion/350314, and notebooks:\nhttps://www.kaggle.com/code/alexandervc/mmscel-cell-cycle-03b-daybydaychange-allcelltypes, https://www.kaggle.com/code/alexandervc/mmscel-cell-cycle-03b-daybydaychange-allcelltypes\n\nAs in the previous scripts we see:\n\n1) Proliferation activity decreases with time\n\n2) Different cell types shows slightly different patterns of proliferation.\n\n3) Donors do not create big differences \n\n#### Technicalities:\nData is big - (105942,  23418) ,  we use sparse fmt created by SBUNZINI https://www.kaggle.com/datasets/sbunzini/open-problems-msci-multiome-sparse-matrices \n\nMetadata. Dataframe with metadata for mutliome train part is prepared. (Using h5py loading packages from .h5 files).\n\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-09-16T19:35:21.784045Z","iopub.execute_input":"2022-09-16T19:35:21.784537Z","iopub.status.idle":"2022-09-16T19:35:21.804879Z","shell.execute_reply.started":"2022-09-16T19:35:21.784500Z","shell.execute_reply":"2022-09-16T19:35:21.803811Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!pip install scanpy\nimport scanpy as sc\nimport anndata\n\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\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\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-09-16T19:35:21.806494Z","iopub.execute_input":"2022-09-16T19:35:21.806835Z","iopub.status.idle":"2022-09-16T19:35:51.912595Z","shell.execute_reply.started":"2022-09-16T19:35:21.806807Z","shell.execute_reply":"2022-09-16T19:35:51.911233Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Create metadata dataframe - cell types, donors, etc","metadata":{}},{"cell_type":"code","source":"%%time\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\nprint()\nprint('Look on TARGETs: ')\nfilename = '/kaggle/input/open-problems-multimodal/train_multi_targets.h5'\nf = h5py.File(filename,'r')#, mode)\nprint(f.keys() )\n\nbarcodes = [t.decode() for t in f['train_multi_targets']['axis1'] ]\nprint(barcodes[:15])\nprint(barcodes[-15:])\n\ndf_meta = pd.DataFrame(index = barcodes)\ndf_meta['Index in Multiome'] = range(len(df_meta))\nprint(df_meta.shape)\n#obs = pd.DataFrame(index = barcodes)\nd = df_cell.reset_index().set_index('cell_id') \nd.columns = ['Index in df_cell'] + list(d.columns[1:])\ndf_meta = df_meta.join(d, how = 'left', )\ndf_meta['str_donor'] = df_meta['donor'].apply(lambda x: str(x))  # will be useful for visualizations\nprint(df_meta.shape)\ndisplay(df_meta)","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:35:51.915319Z","iopub.execute_input":"2022-09-16T19:35:51.915757Z","iopub.status.idle":"2022-09-16T19:36:10.967808Z","shell.execute_reply.started":"2022-09-16T19:35:51.915705Z","shell.execute_reply":"2022-09-16T19:36:10.966622Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Load data","metadata":{}},{"cell_type":"code","source":"%%time\nimport scipy.sparse \n#fn = '/kaggle/input/open-problems-msci-multiome-sparse-matrices/train_multiome_input_sparse.npz'\nfn = '/kaggle/input/open-problems-msci-multiome-sparse-matrices/train_multi_targets_sparse.npz'\n\nY = scipy.sparse.load_npz(fn)# \"../input/multimodal-single-cell-as-sparse-matrix/train_multi_inputs_values.sparse.npz\")\nprint(type(Y) , Y.shape )","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:36:10.969617Z","iopub.execute_input":"2022-09-16T19:36:10.970512Z","iopub.status.idle":"2022-09-16T19:36:27.158349Z","shell.execute_reply.started":"2022-09-16T19:36:10.970477Z","shell.execute_reply":"2022-09-16T19:36:27.157265Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('4.7G used')","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:36:27.160700Z","iopub.execute_input":"2022-09-16T19:36:27.161032Z","iopub.status.idle":"2022-09-16T19:36:27.166180Z","shell.execute_reply.started":"2022-09-16T19:36:27.161003Z","shell.execute_reply":"2022-09-16T19:36:27.165190Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"Y.shape","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:37:28.390642Z","iopub.execute_input":"2022-09-16T19:37:28.391066Z","iopub.status.idle":"2022-09-16T19:37:28.397525Z","shell.execute_reply.started":"2022-09-16T19:37:28.391034Z","shell.execute_reply":"2022-09-16T19:37:28.396834Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Sanity check - sparse data same as original ","metadata":{"execution":{"iopub.status.busy":"2022-09-16T17:39:21.365314Z","iopub.execute_input":"2022-09-16T17:39:21.365760Z","iopub.status.idle":"2022-09-16T17:39:21.371832Z","shell.execute_reply.started":"2022-09-16T17:39:21.365727Z","shell.execute_reply":"2022-09-16T17:39:21.371054Z"}}},{"cell_type":"code","source":"%%time\nSTART = int(0)\nSTOP = START+7# 4_000 # Bigger chunk will crash 16G memory at plt.imshow()\n\ndf_multi_train_y = pd.read_hdf(FP_MULTIOME_TRAIN_TARGETS,start=START,stop=STOP)\ndisplay( df_multi_train_y.head() )\nprint(list( df_multi_train_y.sum(axis = 1)) )\nprint( np.asarray(Y.sum(axis = 1 )).ravel() [:7])\n\nSTART = int(105942+1-7)\nSTOP = START+7# 4_000 # Bigger chunk will crash 16G memory at plt.imshow()\n\ndf_multi_train_y = pd.read_hdf(FP_MULTIOME_TRAIN_TARGETS,start=START,stop=STOP)\ndisplay( df_multi_train_y.head() )\nprint(list( np.round(df_multi_train_y.sum(axis = 1), 3) ))\nprint( np.asarray(Y.sum(axis = 1 )).ravel() [-6:])","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:36:27.167737Z","iopub.execute_input":"2022-09-16T19:36:27.168135Z","iopub.status.idle":"2022-09-16T19:36:27.811444Z","shell.execute_reply.started":"2022-09-16T19:36:27.168103Z","shell.execute_reply":"2022-09-16T19:36:27.810355Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"list_genes = list(df_multi_train_y.columns)\nprint(list_genes[:10])\nprint(list_genes[-10:])","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:36:27.812915Z","iopub.execute_input":"2022-09-16T19:36:27.813220Z","iopub.status.idle":"2022-09-16T19:36:27.821645Z","shell.execute_reply.started":"2022-09-16T19:36:27.813192Z","shell.execute_reply":"2022-09-16T19:36:27.820183Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Update meta data","metadata":{}},{"cell_type":"code","source":"df_meta['n_nonzeros'] = Y.indptr[1:] - Y.indptr[:-1]\ndf_meta\n","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:36:27.823082Z","iopub.execute_input":"2022-09-16T19:36:27.823573Z","iopub.status.idle":"2022-09-16T19:36:27.848319Z","shell.execute_reply.started":"2022-09-16T19:36:27.823532Z","shell.execute_reply":"2022-09-16T19:36:27.847569Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Prepare for cell cycle analysis","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'}\n\n\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\n\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\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","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:36:27.850205Z","iopub.execute_input":"2022-09-16T19:36:27.850596Z","iopub.status.idle":"2022-09-16T19:36:27.870115Z","shell.execute_reply.started":"2022-09-16T19:36:27.850560Z","shell.execute_reply":"2022-09-16T19:36:27.869244Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"s = set(list_genes) & set(G1S_genes_Tirosh)\nprint( len(s), len(set(G1S_genes_Tirosh)) )\ns = set(list_genes) & set(G2M_genes_Tirosh)\nprint(len(s), len(set(G2M_genes_Tirosh)) )\n\n","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:36:27.871198Z","iopub.execute_input":"2022-09-16T19:36:27.871635Z","iopub.status.idle":"2022-09-16T19:36:27.889417Z","shell.execute_reply.started":"2022-09-16T19:36:27.871594Z","shell.execute_reply":"2022-09-16T19:36:27.888624Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# G1/S-G2/M plot for all cells","metadata":{}},{"cell_type":"code","source":"ll = G1S_genes_Tirosh\nI = np.where(pd.Series(list_genes).isin(ll) > 0 )[0] # .sum()\nprint(len(I))\nv1 = np.asarray(Y[:,I].mean(axis = 1)).ravel()\n\nll = G2M_genes_Tirosh\nI = np.where(pd.Series(list_genes).isin(ll) > 0 )[0] # .sum()\nprint(len(I))\nv2 = np.asarray(Y[:,I].mean(axis = 1)).ravel()\n\nxlabel = 'G1/S score'\nylabel = 'G2/M score'\n\nplt.figure(figsize = (20,12))\nmask = np.ones(len(v1)).astype(bool)\nax = sns.scatterplot(x = v1[mask], y = v2[mask], hue = df_meta['n_nonzeros'][mask], palette = 'rainbow' )\nplt.setp(ax.get_legend().get_texts(), fontsize = 20) # for legend text\nplt.setp(ax.get_legend().get_title(), fontsize = 20) # for legend title        \nplt.title('Multiome targets (scRNA-seq). n_cells='+str(mask.sum()),fontsize = 20)\n\nplt.xlabel(xlabel , fontsize = 20)\nplt.ylabel(ylabel, fontsize = 20 )\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:36:27.892128Z","iopub.execute_input":"2022-09-16T19:36:27.892471Z","iopub.status.idle":"2022-09-16T19:36:33.121507Z","shell.execute_reply.started":"2022-09-16T19:36:27.892444Z","shell.execute_reply":"2022-09-16T19:36:33.120133Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"xlabel = 'G1/S score'\nylabel = 'G2/M score'\n\nplt.figure(figsize = (20,12))\nmask = np.ones(len(v1)).astype(bool)\nax = sns.scatterplot(x = v1[mask], y = v2[mask], hue = df_meta['day'][mask], palette = 'rainbow' )\nplt.setp(ax.get_legend().get_texts(), fontsize = 20) # for legend text\nplt.setp(ax.get_legend().get_title(), fontsize = 20) # for legend title    \nplt.title('Multiome targets (scRNA-seq). n_cells='+str(mask.sum()),fontsize = 20)\nplt.xlabel(xlabel , fontsize = 20)\nplt.ylabel(ylabel, fontsize = 20 )\nplt.show()\n\nplt.figure(figsize = (20,12))\nmask = np.ones(len(v1)).astype(bool)\nax = sns.scatterplot(x = v1[mask], y = v2[mask], hue = df_meta['cell_type'][mask], palette = 'rainbow' )\nplt.setp(ax.get_legend().get_texts(), fontsize = 20) # for legend text\nplt.setp(ax.get_legend().get_title(), fontsize = 20) # for legend title     \nplt.title('Multiome targets (scRNA-seq). n_cells='+str(mask.sum()),fontsize = 20)\nplt.xlabel(xlabel , fontsize = 20)\nplt.ylabel(ylabel, fontsize = 20 )\nplt.show()\n\nplt.figure(figsize = (20,12))\nmask = np.ones(len(v1)).astype(bool)\nax = sns.scatterplot(x = v1[mask], y = v2[mask], hue = df_meta['str_donor'][mask], palette = 'rainbow' , alpha = 0.1)\nplt.setp(ax.get_legend().get_texts(), fontsize = 20) # for legend text\nplt.setp(ax.get_legend().get_title(), fontsize = 20) # for legend title     \nplt.title('Multiome targets (scRNA-seq). n_cells='+str(mask.sum()),fontsize = 20)\nplt.xlabel(xlabel , fontsize = 20)\nplt.ylabel(ylabel, fontsize = 20 )\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:36:33.123049Z","iopub.execute_input":"2022-09-16T19:36:33.123340Z","iopub.status.idle":"2022-09-16T19:36:42.751122Z","shell.execute_reply.started":"2022-09-16T19:36:33.123315Z","shell.execute_reply":"2022-09-16T19:36:42.750108Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print( df_meta['day'].unique() )\nprint( df_meta['cell_type'].unique() )\ndf_meta['cell_type'].value_counts()","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:36:42.752309Z","iopub.execute_input":"2022-09-16T19:36:42.752630Z","iopub.status.idle":"2022-09-16T19:36:42.774691Z","shell.execute_reply.started":"2022-09-16T19:36:42.752601Z","shell.execute_reply":"2022-09-16T19:36:42.773680Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# G1/S-G2/M plot loop over cell types and days","metadata":{}},{"cell_type":"code","source":"\nfor uv1 in ['NeuP', 'HSC', 'MasP', 'MkP', 'EryP', 'MoP', 'BP']:\n    col = 'cell_type'\n\n    mask1 = (df_meta[col] == uv1)\n    str_inf1 = 'Multiome targets  (scRNA-seq). ' + col + '= ' +str(uv1) + ' ' + str(dict_cell_types[uv1])\n    print()\n    print(uv1)\n    print(str_inf1, mask1.sum() )\n    \n    for uv in [2, 3, 4, 7]:\n        #uv = 2\n        col = 'day'\n        mask = mask1 & (df_meta[col] == uv)\n        str_inf = str_inf1 + '  ' + col + '= ' +str(uv)  + ' n_cells= ' + str(mask.sum() )\n        print(str_inf, mask.sum() )\n        print()\n\n        plt.figure(figsize = (20,12))\n        ax = sns.scatterplot(x = v1[mask], y = v2[mask], hue = df_meta['n_nonzeros'][mask], palette = 'rainbow' )\n        plt.setp(ax.get_legend().get_texts(), fontsize = 20) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize = 20) # for legend title        \n        plt.title(str_inf,fontsize = 20)\n        plt.xlabel(xlabel , fontsize = 20)\n        plt.ylabel(ylabel, fontsize = 20 )\n        plt.show()\n","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:36:42.775720Z","iopub.execute_input":"2022-09-16T19:36:42.775991Z","iopub.status.idle":"2022-09-16T19:36:57.021412Z","shell.execute_reply.started":"2022-09-16T19:36:42.775966Z","shell.execute_reply":"2022-09-16T19:36:57.020449Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# G1/S-G2/M plot for day 7, loop over cell types ","metadata":{}},{"cell_type":"code","source":"\nfor uv1 in [ 7]:\n    #uv = 2\n    col = 'day'\n    mask1 = (df_meta[col] == uv1)\n    str_inf1 = 'Multiome targets  (scRNA-seq). ' + col + '= ' +str(uv1)\n    print()\n    print(str_inf1, mask1.sum() )\n    \n    col = 'cell_type'\n    for uv in ['NeuP', 'HSC', 'MasP', 'MkP', 'EryP', 'MoP', 'BP']:\n        mask = mask1 & (df_meta[col] == uv)\n        str_inf = str_inf1 + '  ' + col + '= ' +str(uv)+ ' ' + str(dict_cell_types[uv]) + ' n_cells= ' + str(mask.sum() )\n        print(str_inf, mask.sum() )\n\n        plt.figure(figsize = (20,12))\n        ax = sns.scatterplot(x = v1[mask], y = v2[mask], hue = df_meta['n_nonzeros'][mask], palette = 'rainbow' )\n        plt.setp(ax.get_legend().get_texts(), fontsize = 20) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize = 20) # for legend title        \n        plt.title(str_inf,fontsize = 20)\n        plt.xlabel(xlabel , fontsize = 20)\n        plt.ylabel(ylabel, fontsize = 20 )\n        plt.show()\n","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:36:57.022480Z","iopub.execute_input":"2022-09-16T19:36:57.022784Z","iopub.status.idle":"2022-09-16T19:37:00.336234Z","shell.execute_reply.started":"2022-09-16T19:36:57.022757Z","shell.execute_reply":"2022-09-16T19:37:00.335268Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# UMAP for cell cycle - not quite succesful as probably should be expected. \n\nThe cyclic structure does not magically appear, thus umap since is not acting as denoising (MAGIC etc).\n","metadata":{}},{"cell_type":"code","source":"import umap\nreducer = umap.UMAP()\n","metadata":{"execution":{"iopub.status.busy":"2022-09-16T20:25:24.203292Z","iopub.execute_input":"2022-09-16T20:25:24.203712Z","iopub.status.idle":"2022-09-16T20:25:41.951274Z","shell.execute_reply.started":"2022-09-16T20:25:24.203679Z","shell.execute_reply":"2022-09-16T20:25:41.950087Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\nfor uv1 in [ 7]:\n    #uv = 2\n    col = 'day'\n    mask1 = (df_meta[col] == uv1)\n    str_inf1 = 'Multiome targets  (scRNA-seq). ' + col + '= ' +str(uv1)\n    print()\n    print(str_inf1, mask1.sum() )\n    \n    col = 'cell_type'\n    for uv in ['NeuP', 'HSC', 'MasP', 'MkP', 'EryP', 'MoP', 'BP']:\n        mask = mask1 & (df_meta[col] == uv)\n        str_inf = str_inf1 + '  ' + col + '= ' +str(uv)+ ' ' + str(dict_cell_types[uv]) + ' n_cells= ' + str(mask.sum() )\n        print(str_inf, mask.sum() )\n\n        I = np.where(pd.Series(list_genes).isin(G1S_genes_Tirosh + G2M_genes_Tirosh) > 0 )[0] # .sum()\n        print(len(I))\n        X = Y[mask][:,I]#  Y[mask,I ]# np.asarray(Y[:,I].mean(axis = 1)).ravel()\n\n        r = reducer.fit_transform(X)\n\n        plt.figure(figsize = (20,12))\n        #ax = sns.scatterplot(x = v1[mask], y = v2[mask], hue = df_meta['n_nonzeros'][mask], palette = 'rainbow' )\n        ax = sns.scatterplot(x = r[:,0], y = r[:,1], hue = df_meta['n_nonzeros'][mask], palette = 'rainbow' )\n        plt.setp(ax.get_legend().get_texts(), fontsize = 20) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize = 20) # for legend title        \n        plt.title(str_inf,fontsize = 20)\n        plt.xlabel('UMAP1' , fontsize = 20)\n        plt.ylabel('UMAP2', fontsize = 20 )\n        plt.show()\n","metadata":{"execution":{"iopub.status.busy":"2022-09-16T20:30:47.963156Z","iopub.execute_input":"2022-09-16T20:30:47.963776Z","iopub.status.idle":"2022-09-16T20:32:10.750078Z","shell.execute_reply.started":"2022-09-16T20:30:47.963707Z","shell.execute_reply":"2022-09-16T20:32:10.748856Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Changing params of umap ","metadata":{}},{"cell_type":"code","source":"for n_neighbors in (5, 15, 50, 100): # min_dist=0.1\n    for min_dist in [0.1,0.5, 0.9]:\n        reducer = umap.UMAP(n_neighbors=n_neighbors,  min_dist=min_dist )\n        for uv1 in [ 7]:\n            #uv = 2\n            col = 'day'\n            mask1 = (df_meta[col] == uv1)\n            str_inf1 ='UMAP' + ' n_neighbors='+str(n_neighbors) +  ' min_dist='+str(min_dist) + ' Multiome targets  (scRNA-seq). ' + col + '= ' +str(uv1)\n            print()\n            print(str_inf1, mask1.sum() )\n\n            col = 'cell_type'\n            for uv in [ 'MasP',]:# ['NeuP', 'HSC', 'MasP', 'MkP', 'EryP', 'MoP', 'BP']:\n                mask = mask1 & (df_meta[col] == uv)\n                str_inf = str_inf1 + '  ' + col + '= ' +str(uv)+ ' ' + str(dict_cell_types[uv]) + ' n_cells= ' + str(mask.sum() )\n                print(str_inf, mask.sum() )\n\n                I = np.where(pd.Series(list_genes).isin(G1S_genes_Tirosh + G2M_genes_Tirosh) > 0 )[0] # .sum()\n                print(len(I))\n                X = Y[mask][:,I]#  Y[mask,I ]# np.asarray(Y[:,I].mean(axis = 1)).ravel()\n\n                r = reducer.fit_transform(X)\n\n                plt.figure(figsize = (20,12))\n                #ax = sns.scatterplot(x = v1[mask], y = v2[mask], hue = df_meta['n_nonzeros'][mask], palette = 'rainbow' )\n                ax = sns.scatterplot(x = r[:,0], y = r[:,1], hue = df_meta['n_nonzeros'][mask], palette = 'rainbow' )\n                plt.setp(ax.get_legend().get_texts(), fontsize = 20) # for legend text\n                plt.setp(ax.get_legend().get_title(), fontsize = 20) # for legend title        \n                plt.title(str_inf,fontsize = 20)\n                plt.xlabel('UMAP1' , fontsize = 20)\n                plt.ylabel('UMAP2', fontsize = 20 )\n                plt.show()\n","metadata":{"execution":{"iopub.status.busy":"2022-09-16T20:42:03.352787Z","iopub.execute_input":"2022-09-16T20:42:03.353231Z","iopub.status.idle":"2022-09-16T20:47:28.639263Z","shell.execute_reply.started":"2022-09-16T20:42:03.353200Z","shell.execute_reply":"2022-09-16T20:47:28.638112Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('%.1f seconds passed total '%(time.time()-t0start) )\n","metadata":{"execution":{"iopub.status.busy":"2022-09-16T19:37:00.337621Z","iopub.execute_input":"2022-09-16T19:37:00.338008Z","iopub.status.idle":"2022-09-16T19:37:00.344012Z","shell.execute_reply.started":"2022-09-16T19:37:00.337970Z","shell.execute_reply":"2022-09-16T19:37:00.342865Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{"execution":{"iopub.status.busy":"2022-09-16T20:30:34.348253Z","iopub.execute_input":"2022-09-16T20:30:34.349212Z","iopub.status.idle":"2022-09-16T20:30:34.522352Z","shell.execute_reply.started":"2022-09-16T20:30:34.349173Z","shell.execute_reply":"2022-09-16T20:30:34.521207Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}