{"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\nWe consider count non-zero elements in ATAC-seq (features of Multiome task); sum of ATAC-seq values and PCAs from ATAC seq.\n\n**Findings Briefly:** \n\n    1) day dependence for sums and counts\n            counts of non-zeroes sharply grow at day 3, and little decrease at day 7\n            sum decrease at day 3, and little increase at day 2\n            so day 3 is quite distinguished , in contrast to day 7 for scRNA-seq data\n    2) PCA8 from ATAC is weakly correlated to GEX features like sums and cell cycle\n    3) Sums of ATAC values and PCA1 are weakly related to cell proliferation activity (calculated by GEX side) see visualizations below \n\n**Remark:** see also discussion: https://www.kaggle.com/competitions/open-problems-multimodal/discussion/353629\nand notebook: https://www.kaggle.com/alexandervc/mmscel-count-nonzero-genes-decrease-daily (similar analysis for scRNA-seq - Multiome targets and CITE-seq features). \n\n**Previous EDAs** on ATAC-seq part of Multiome: \n\n1) https://www.kaggle.com/code/alexandervc/mmscel-atac-seq-eda-01-dimred-visualizations\n    \n    Conclusions: \n    a) cell types can be seen - quite well. (PCA2,3)\n    b) PCA1 - strongly correlated with total sum - probably just technical (at least partly) variation.\n    \n2) https://www.kaggle.com/code/alexandervc/heatmap-of-atac-seq-multiome-data\n\n    Conclusions:  (Heatmaps created and plt.spy used for heatmaps for sparse matrices) \n    a) Main patterns - zero or not - are not quite changing from cell to cells - we see vertical white / blue strips\n    b) Day patterns can be seen just on heatmaps (plt.spy - special heatmaps for sparse matrices)\n\n\n**Remark2:** 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. See also discussion https://www.kaggle.com/competitions/open-problems-multimodal/discussion/350314, and notebooks: https://www.kaggle.com/code/alexandervc/mmscel-cell-cycle-03b-daybydaychange-allcelltypes, https://www.kaggle.com/code/alexandervc/mmscel-cell-cycle-03b-daybydaychange-allcelltypes\nhttps://www.kaggle.com/code/alexandervc/mmscel-eda-multiome-targets-01-cell-cycle","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-19T18:35:40.857634Z","iopub.execute_input":"2022-09-19T18:35:40.858168Z","iopub.status.idle":"2022-09-19T18:35:40.900621Z","shell.execute_reply.started":"2022-09-19T18:35:40.858063Z","shell.execute_reply":"2022-09-19T18:35:40.898968Z"},"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\n# df_cell = pd.read_csv(FP_CELL_METADATA)\n# df_cell","metadata":{"execution":{"iopub.status.busy":"2022-09-19T18:35:41.564455Z","iopub.execute_input":"2022-09-19T18:35:41.565164Z","iopub.status.idle":"2022-09-19T18:36:25.982725Z","shell.execute_reply.started":"2022-09-19T18:35:41.565122Z","shell.execute_reply":"2022-09-19T18:36:25.980890Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Load ATAC-seq data in sparse format (created by SBUNZINI )","metadata":{}},{"cell_type":"code","source":"%%time\nstr_data_inf = ' ATAC-seq - Train part of Multiome '\nimport scipy.sparse \nfn = '/kaggle/input/open-problems-msci-multiome-sparse-matrices/train_multiome_input_sparse.npz'\ntrain_inputs = scipy.sparse.load_npz(fn)# \"../input/multimodal-single-cell-as-sparse-matrix/train_multi_inputs_values.sparse.npz\")\n\nprint('8.7G consumed')","metadata":{"execution":{"iopub.status.busy":"2022-09-19T18:36:25.988316Z","iopub.execute_input":"2022-09-19T18:36:25.989735Z","iopub.status.idle":"2022-09-19T18:37:45.791152Z","shell.execute_reply.started":"2022-09-19T18:36:25.989670Z","shell.execute_reply":"2022-09-19T18:37:45.789698Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_inputs","metadata":{"execution":{"iopub.status.busy":"2022-09-19T18:37:45.792805Z","iopub.execute_input":"2022-09-19T18:37:45.793207Z","iopub.status.idle":"2022-09-19T18:37:45.803790Z","shell.execute_reply.started":"2022-09-19T18:37:45.793168Z","shell.execute_reply":"2022-09-19T18:37:45.802534Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('8.7G RAM consumed')","metadata":{"execution":{"iopub.status.busy":"2022-09-19T18:37:45.806173Z","iopub.execute_input":"2022-09-19T18:37:45.806537Z","iopub.status.idle":"2022-09-19T18:37:45.816725Z","shell.execute_reply.started":"2022-09-19T18:37:45.806501Z","shell.execute_reply":"2022-09-19T18:37:45.815537Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Load meta data and scores, update it with ATAC engeneered features","metadata":{}},{"cell_type":"code","source":"f = '/kaggle/input/scrnaseq-competition-open-problems-adata-fmt/multiome_meta_and_scores_etc.csv'\ndf_meta = pd.read_csv(f, index_col = 0)\ndf_meta","metadata":{"execution":{"iopub.status.busy":"2022-09-19T18:39:24.723040Z","iopub.execute_input":"2022-09-19T18:39:24.723577Z","iopub.status.idle":"2022-09-19T18:39:25.830575Z","shell.execute_reply.started":"2022-09-19T18:39:24.723535Z","shell.execute_reply":"2022-09-19T18:39:25.829262Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nv = train_inputs.sum(axis = 1)\ndf_meta['sum ATAC'] = np.asarray(v).ravel()\n","metadata":{"execution":{"iopub.status.busy":"2022-09-19T18:39:40.850070Z","iopub.execute_input":"2022-09-19T18:39:40.850475Z","iopub.status.idle":"2022-09-19T18:39:41.154429Z","shell.execute_reply.started":"2022-09-19T18:39:40.850439Z","shell.execute_reply":"2022-09-19T18:39:41.153151Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\ndf_meta['n_nonzeros ATAC'] = train_inputs.indptr[1:] - train_inputs.indptr[:-1] # Number of nonzero elements in sparse matrix row by row. Borrowed from SO: \ndf_meta.head(1)","metadata":{"execution":{"iopub.status.busy":"2022-09-19T18:39:42.412079Z","iopub.execute_input":"2022-09-19T18:39:42.412755Z","iopub.status.idle":"2022-09-19T18:39:42.442106Z","shell.execute_reply.started":"2022-09-19T18:39:42.412714Z","shell.execute_reply":"2022-09-19T18:39:42.440927Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Look on count non-zero ATAC elements - see clear dependence on day with day 3 significantly larger ","metadata":{}},{"cell_type":"code","source":"plt.figure(figsize = (20,5))\n#df_meta['n_nonzeros'].plot()\nplt.plot( df_meta['n_nonzeros ATAC'].values, label = 'n_nonzeros ATAC' )\nv = df_meta['day']\nplt.plot(v.values.ravel()*1000, label = 'Day')\nplt.plot( df_meta['donor'].values/1, label = 'donor')\nplt.legend(fontsize = 14)\nplt.suptitle(str_data_inf , fontsize = 20)\nplt.title('Count of nonzero ATAC elements per cell shows clear daily dependent pattern. Day 3 (second probe day) is significantly larger ')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-09-19T18:39:47.153154Z","iopub.execute_input":"2022-09-19T18:39:47.153588Z","iopub.status.idle":"2022-09-19T18:39:49.236086Z","shell.execute_reply.started":"2022-09-19T18:39:47.153551Z","shell.execute_reply":"2022-09-19T18:39:49.234759Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"d = np.round(df_meta.groupby(['donor','day'])['n_nonzeros ATAC','sum ATAC','sum GEX', 'n_nonzeros GEX','G1S+G2M cell cycle score'].mean())\ndisplay(d)\nd2 = np.round(df_meta.groupby(['day'])['n_nonzeros ATAC','sum ATAC','sum GEX', 'n_nonzeros GEX','G1S+G2M cell cycle score'].mean())\ndisplay(d2)\nplt.bar(x = d2.index , height= d2['sum ATAC'].values)#.plot()","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:37:06.279237Z","iopub.execute_input":"2022-09-19T19:37:06.279691Z","iopub.status.idle":"2022-09-19T19:37:06.546769Z","shell.execute_reply.started":"2022-09-19T19:37:06.279652Z","shell.execute_reply":"2022-09-19T19:37:06.545504Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"d2 = np.round(df_meta.groupby(['cell_type'])['n_nonzeros ATAC','sum ATAC','sum GEX', 'n_nonzeros GEX','G1S+G2M cell cycle score'].mean())\ndisplay(d2)\nd2 = np.round(df_meta.groupby(['day','cell_type'])['n_nonzeros ATAC','sum ATAC','sum GEX', 'n_nonzeros GEX','G1S+G2M cell cycle score'].mean())\ndisplay(d2)\n","metadata":{"execution":{"iopub.status.busy":"2022-09-19T20:01:57.827874Z","iopub.execute_input":"2022-09-19T20:01:57.828409Z","iopub.status.idle":"2022-09-19T20:01:57.920290Z","shell.execute_reply.started":"2022-09-19T20:01:57.828363Z","shell.execute_reply":"2022-09-19T20:01:57.919046Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Look on ATAC seq sum of values - day 3 shows decreas ","metadata":{}},{"cell_type":"code","source":"plt.figure(figsize = (20,5))\n#df_meta['n_nonzeros'].plot()\nplt.plot( df_meta['sum ATAC'].values, label = 'sum ATAC' )\nv = df_meta['day']\nplt.plot(v.values.ravel()*1000, label = 'Day')\nplt.plot( df_meta['donor'].values/10, label = 'donor')\nplt.legend(fontsize = 14)\nplt.suptitle(str_data_inf , fontsize = 20)\nplt.title('Sum of  ATAC values  - shows  daily dependent pattern. Day 3 (second probe day) is smaller ')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-09-19T17:56:32.238959Z","iopub.execute_input":"2022-09-19T17:56:32.239455Z","iopub.status.idle":"2022-09-19T17:56:34.555190Z","shell.execute_reply.started":"2022-09-19T17:56:32.239413Z","shell.execute_reply":"2022-09-19T17:56:34.553967Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Truncated SVD dimensional reduction (one donor - to get result faster)","metadata":{}},{"cell_type":"code","source":"from sklearn.decomposition import TruncatedSVD","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:07:09.561738Z","iopub.execute_input":"2022-09-19T19:07:09.562213Z","iopub.status.idle":"2022-09-19T19:07:09.760413Z","shell.execute_reply.started":"2022-09-19T19:07:09.562169Z","shell.execute_reply":"2022-09-19T19:07:09.759174Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nfor donor in [32606]:# df_meta_multiome_train['donor'].unique():\n    #donor = 32606\n    #day = 2\n    str_inf = ''\n    str_inf = ' donor= ' + str(donor) \n    #str_inf += ' day= ' + str(day)\n    #mask = (df_meta_multiome_train['day'] == day) # & (df_meta_multiome_train['donor'] == donor)\n    mask = (df_meta['donor'] == donor)\n\n    print(mask.sum(), str_inf )\n    tsvd = TruncatedSVD(n_components = 32,random_state = 4223)\n    t0 = time.time()\n    r = tsvd.fit_transform(train_inputs[mask,:])\n    print('%.1f'%(time.time()-t0) + ' seconds passed. TruncatedSVD finished.')\n","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:07:10.106816Z","iopub.execute_input":"2022-09-19T19:07:10.107256Z","iopub.status.idle":"2022-09-19T19:10:11.948922Z","shell.execute_reply.started":"2022-09-19T19:07:10.107217Z","shell.execute_reply":"2022-09-19T19:10:11.947686Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"tsvd.components_.shape","metadata":{"execution":{"iopub.status.busy":"2022-09-19T20:48:02.609974Z","iopub.execute_input":"2022-09-19T20:48:02.610653Z","iopub.status.idle":"2022-09-19T20:48:02.627051Z","shell.execute_reply.started":"2022-09-19T20:48:02.610600Z","shell.execute_reply":"2022-09-19T20:48:02.621087Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for i in range(12):\n    plt.figure(figsize = (20,4))\n    plt.plot(tsvd.components_[i,:1000]) # .shape\n    plt.title('PCA '+str(i+1) + ' 1000 elements')\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2022-09-19T20:52:43.046357Z","iopub.execute_input":"2022-09-19T20:52:43.046793Z","iopub.status.idle":"2022-09-19T20:52:45.807479Z","shell.execute_reply.started":"2022-09-19T20:52:43.046757Z","shell.execute_reply":"2022-09-19T20:52:45.806097Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Analysis of correlations for RNAseq (GEX)  and ATAC engineered features. PCA8 from ATAC shows some correlations with RNAseq features ","metadata":{"execution":{"iopub.status.busy":"2022-09-19T18:11:04.513449Z","iopub.execute_input":"2022-09-19T18:11:04.514016Z","iopub.status.idle":"2022-09-19T18:11:04.520939Z","shell.execute_reply.started":"2022-09-19T18:11:04.513973Z","shell.execute_reply":"2022-09-19T18:11:04.519565Z"}}},{"cell_type":"code","source":"for i in range(r.shape[1]):\n    df_meta.loc[mask, 'ATAC PCA'+str(i+1)] =r[:,i]\ncm = df_meta[mask].corr()\n","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:10:55.096921Z","iopub.execute_input":"2022-09-19T19:10:55.097442Z","iopub.status.idle":"2022-09-19T19:10:55.612746Z","shell.execute_reply.started":"2022-09-19T19:10:55.097397Z","shell.execute_reply":"2022-09-19T19:10:55.611410Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize = (20,14))\nsns.heatmap(cm)","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:10:58.602258Z","iopub.execute_input":"2022-09-19T19:10:58.602767Z","iopub.status.idle":"2022-09-19T19:11:01.219963Z","shell.execute_reply.started":"2022-09-19T19:10:58.602721Z","shell.execute_reply":"2022-09-19T19:11:01.218955Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import pandas as pd\npd.set_option('display.max_rows', 500)\npd.set_option('display.max_columns', 500)\npd.set_option('display.width', 1000)\n\ncm=np.round(cm,2)\ncm","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:11:03.044311Z","iopub.execute_input":"2022-09-19T19:11:03.044740Z","iopub.status.idle":"2022-09-19T19:11:03.402947Z","shell.execute_reply.started":"2022-09-19T19:11:03.044703Z","shell.execute_reply":"2022-09-19T19:11:03.401749Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# PCA8 PCA* plots showing some dependence between GEX features and ATAC PCA8 ","metadata":{}},{"cell_type":"code","source":"for (i,j) in [ (0,7),(1,7), (2,7),(3,7),(4,7), (5,7),(6,7),(8,7),(9,7)]:\n    for col in ['G1S+G2M cell cycle score']:# ,  'sum GEX']:\n        plt.figure(figsize = (20,10))\n        ax = sns.scatterplot(x = r[:,i], y = r[:,j], hue = df_meta[col][mask], palette = 'rainbow', alpha = 0.9)\n        if 1:\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('ATAC-seq. Train. PCA (TruncSVD). ' + str_inf  +   ' n_cells='+str(len(r)), fontsize = 20)\n        plt.xlabel('PCA'+str(i+1), fontsize = 20)\n        plt.ylabel('PCA'+str(j+1), fontsize = 20)\n        plt.show() ","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:11:20.160765Z","iopub.execute_input":"2022-09-19T19:11:20.161282Z","iopub.status.idle":"2022-09-19T19:11:32.998011Z","shell.execute_reply.started":"2022-09-19T19:11:20.161240Z","shell.execute_reply":"2022-09-19T19:11:32.996908Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# GEX features (sums, cell cycle) are weakly related to ATAC PCA1 (similar to sum)","metadata":{}},{"cell_type":"code","source":"for (i,j) in [ (0,1)]: #,(0,2),(1,2), (0,3),(1,3),(2,3), (0,4),(1,4),(2,4),(3,4)]:\n    for col in ['n_nonzeros ATAC', 'sum ATAC', 'n_nonzeros GEX', #'G1S score', 'G2M score',\n       'G1S+G2M cell cycle score',  'FAST score', 'sum GEX']:# , 'Sum chrY', 'Sum chrMT']:\n        plt.figure(figsize = (20,10))\n        ax = sns.scatterplot(x = r[:,i], y = r[:,j], hue = df_meta[col][mask], palette = 'rainbow', alpha = 0.9)\n        if 1:\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('ATAC-seq. Train. PCA (TruncSVD). ' + str_inf  +   ' n_cells='+str(len(r)), fontsize = 20)\n        plt.xlabel('PCA'+str(i+1), fontsize = 20)\n        plt.ylabel('PCA'+str(j+1), fontsize = 20)\n        plt.show() ","metadata":{"execution":{"iopub.status.busy":"2022-09-19T20:36:09.372843Z","iopub.execute_input":"2022-09-19T20:36:09.373407Z","iopub.status.idle":"2022-09-19T20:36:18.559569Z","shell.execute_reply.started":"2022-09-19T20:36:09.373354Z","shell.execute_reply":"2022-09-19T20:36:18.557908Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# PCA3,5 where cell types are clearly seen - colored by different GEX and ATAC features - not much is seen","metadata":{}},{"cell_type":"code","source":"for (i,j) in [ (2,4)]: #,(0,2),(1,2), (0,3),(1,3),(2,3), (0,4),(1,4),(2,4),(3,4)]:\n    for col in ['cell_type', 'n_nonzeros GEX', #'G1S score', 'G2M score',\n       'G1S+G2M cell cycle score',  'FAST score', 'sum GEX', 'n_nonzeros ATAC','sum ATAC']:# , 'Sum chrY', 'Sum chrMT']:\n        plt.figure(figsize = (20,10))\n        ax = sns.scatterplot(x = r[:,i], y = r[:,j], hue = df_meta[col][mask], palette = 'rainbow', alpha = 0.9)\n        if 1:\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('ATAC-seq. Train. PCA (TruncSVD). ' + str_inf  +   ' n_cells='+str(len(r)), fontsize = 20)\n        plt.xlabel('PCA'+str(i+1), fontsize = 20)\n        plt.ylabel('PCA'+str(j+1), fontsize = 20)\n        plt.show() ","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:25:44.108482Z","iopub.execute_input":"2022-09-19T19:25:44.108981Z","iopub.status.idle":"2022-09-19T19:25:53.786222Z","shell.execute_reply.started":"2022-09-19T19:25:44.108923Z","shell.execute_reply":"2022-09-19T19:25:53.785039Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# UMAP plots from TruncatedSVD 32 ","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:18:46.616442Z","iopub.execute_input":"2022-09-19T19:18:46.616966Z","iopub.status.idle":"2022-09-19T19:18:46.627461Z","shell.execute_reply.started":"2022-09-19T19:18:46.616922Z","shell.execute_reply":"2022-09-19T19:18:46.625132Z"}}},{"cell_type":"code","source":"%%time\nprint(r.shape)\nimport umap\nreducer = umap.UMAP()\nr2 = reducer.fit_transform(r)\n\nfor (i,j) in [(0,1)]:# ,(0,2),(1,2), (0,3),(1,3),(2,3), (0,4),(1,4),(2,4),(3,4)]:\n    plt.figure(figsize = (20,10))\n    ax = sns.scatterplot(x=r2[:,i],y=r2[:,j], hue = df_meta['day'][mask], palette = 'rainbow')\n    if 1:\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('ATAC-seq. Train. UMAP after PCA (TruncSVD) 10. ' + str_inf  +   ' n_cells='+str(len(r)), fontsize = 20)\n    plt.xlabel('UMAP'+str(i+1), fontsize = 20)\n    plt.ylabel('UMAP'+str(j+1), fontsize = 20)\n    plt.show() ","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:18:48.526017Z","iopub.execute_input":"2022-09-19T19:18:48.526444Z","iopub.status.idle":"2022-09-19T19:19:59.340114Z","shell.execute_reply.started":"2022-09-19T19:18:48.526409Z","shell.execute_reply":"2022-09-19T19:19:59.338767Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for (i,j) in [(0,1)]:# ,(0,2),(1,2), (0,3),(1,3),(2,3), (0,4),(1,4),(2,4),(3,4)]:\n    for col in ['cell_type', 'n_nonzeros GEX', 'G1S score', 'G2M score',\n       'G1S+G2M cell cycle score',  'FAST score', 'sum GEX', 'sum chrY', 'sum chrMT']:\n        plt.figure(figsize = (20,10))\n        ax = sns.scatterplot(x=r2[:,i],y=r2[:,j], hue = df_meta[col][mask], palette = 'rainbow')\n        if 1:\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('ATAC-seq. Train. UMAP after PCA (TruncSVD) 10. ' + str_inf  +   ' n_cells='+str(len(r)), fontsize = 20)\n        plt.xlabel('UMAP'+str(i+1), fontsize = 20)\n        plt.ylabel('UMAP'+str(j+1), fontsize = 20)\n        plt.show()         ","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:22:25.196443Z","iopub.execute_input":"2022-09-19T19:22:25.196891Z","iopub.status.idle":"2022-09-19T19:22:37.855120Z","shell.execute_reply.started":"2022-09-19T19:22:25.196856Z","shell.execute_reply":"2022-09-19T19:22:37.853728Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Save df_meta to disk","metadata":{}},{"cell_type":"code","source":"print(df_meta.shape)\ndf_meta.head(2)","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:54:52.103066Z","iopub.execute_input":"2022-09-19T19:54:52.103862Z","iopub.status.idle":"2022-09-19T19:54:52.162370Z","shell.execute_reply.started":"2022-09-19T19:54:52.103803Z","shell.execute_reply":"2022-09-19T19:54:52.161118Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_meta.to_csv('multiome_meta_and_scores_etc_RNA_ATAC.csv')","metadata":{"execution":{"iopub.status.busy":"2022-09-19T19:56:11.878134Z","iopub.execute_input":"2022-09-19T19:56:11.878622Z","iopub.status.idle":"2022-09-19T19:56:17.016427Z","shell.execute_reply.started":"2022-09-19T19:56:11.878579Z","shell.execute_reply":"2022-09-19T19:56:17.015108Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"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-09-16T17:00:34.816138Z","iopub.execute_input":"2022-09-16T17:00:34.816547Z","iopub.status.idle":"2022-09-16T17:00:34.824518Z","shell.execute_reply.started":"2022-09-16T17:00:34.816512Z","shell.execute_reply":"2022-09-16T17:00:34.822628Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}