{"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":"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":"2023-09-22T19:57:10.950625Z","iopub.execute_input":"2023-09-22T19:57:10.951073Z","iopub.status.idle":"2023-09-22T19:57:11.517436Z","shell.execute_reply.started":"2023-09-22T19:57:10.951036Z","shell.execute_reply":"2023-09-22T19:57:11.516265Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_path = '/kaggle/input/open-problems-single-cell-perturbations/de_train.parquet'\ndf_train = pd.read_parquet(train_path)\nprint(df_train.shape)\ndf_train","metadata":{"execution":{"iopub.status.busy":"2023-09-22T19:57:11.519539Z","iopub.execute_input":"2023-09-22T19:57:11.520331Z","iopub.status.idle":"2023-09-22T19:57:14.995066Z","shell.execute_reply.started":"2023-09-22T19:57:11.520291Z","shell.execute_reply":"2023-09-22T19:57:14.993308Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Get the list of compounds that have data all six cell types\ncompound_list_train = list((df_train['sm_name'].value_counts()).index[df_train['sm_name'].value_counts() == 6])\n\n# Reomve positive controls\ncompound_list_train.remove('Dabrafenib')\ncompound_list_train.remove('Belinostat')\ndisplay(compound_list_train)","metadata":{"execution":{"iopub.status.busy":"2023-09-22T19:57:14.997757Z","iopub.execute_input":"2023-09-22T19:57:14.998652Z","iopub.status.idle":"2023-09-22T19:57:15.020937Z","shell.execute_reply.started":"2023-09-22T19:57:14.998610Z","shell.execute_reply":"2023-09-22T19:57:15.019498Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# filter by the above compound list\nfiltered_df_train = df_train[df_train['sm_name'].isin(compound_list_train)]\n\n# filter by the target cell types\nfiltered_df_train = filtered_df_train[filtered_df_train['cell_type'].isin(['B cells', 'Myeloid cells'])]\nfiltered_df_train","metadata":{"execution":{"iopub.status.busy":"2023-09-22T19:57:15.024704Z","iopub.execute_input":"2023-09-22T19:57:15.025330Z","iopub.status.idle":"2023-09-22T19:57:15.100895Z","shell.execute_reply.started":"2023-09-22T19:57:15.025284Z","shell.execute_reply":"2023-09-22T19:57:15.099349Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# MRRMSE","metadata":{}},{"cell_type":"markdown","source":"$$\\textrm{MRRMSE} = \\frac{1}{R}\\sum_{i=1}^R\\left(\\frac{1}{n} \\sum_{j=1}^{n} (y_{ij} - \\widehat{y}_{ij})^2\\right)^{1/2}$$","metadata":{}},{"cell_type":"code","source":"((filtered_df_train.iloc[:, 5:] ** 2).mean(axis=1) ** 0.5).mean()","metadata":{"execution":{"iopub.status.busy":"2023-09-22T19:57:15.103087Z","iopub.execute_input":"2023-09-22T19:57:15.103573Z","iopub.status.idle":"2023-09-22T19:57:15.125617Z","shell.execute_reply.started":"2023-09-22T19:57:15.103539Z","shell.execute_reply":"2023-09-22T19:57:15.124276Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def mrrmse(df: pd.DataFrame):\n\n    return (\n            (df.iloc[:, 5:] ** 2                    # (y_ij - yhat_{ij}) ^ 2\n            ).mean(axis=1) ** 0.5                   # take the mean of n genes then sqaure root\n           ).mean()                                    # take the mean by R","metadata":{"execution":{"iopub.status.busy":"2023-09-22T19:57:15.127554Z","iopub.execute_input":"2023-09-22T19:57:15.128016Z","iopub.status.idle":"2023-09-22T19:57:15.136049Z","shell.execute_reply.started":"2023-09-22T19:57:15.127967Z","shell.execute_reply":"2023-09-22T19:57:15.134343Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Entire Train Set","metadata":{}},{"cell_type":"code","source":"mrrmse(filtered_df_train)","metadata":{"execution":{"iopub.status.busy":"2023-09-22T19:57:15.138522Z","iopub.execute_input":"2023-09-22T19:57:15.139389Z","iopub.status.idle":"2023-09-22T19:57:15.173443Z","shell.execute_reply.started":"2023-09-22T19:57:15.139353Z","shell.execute_reply":"2023-09-22T19:57:15.172279Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from random import shuffle\nfrom tqdm.notebook import tqdm\nfrom plotly import express as px\nimport seaborn as sns","metadata":{"execution":{"iopub.status.busy":"2023-09-22T20:23:10.502658Z","iopub.execute_input":"2023-09-22T20:23:10.503420Z","iopub.status.idle":"2023-09-22T20:23:11.361200Z","shell.execute_reply.started":"2023-09-22T20:23:10.503382Z","shell.execute_reply":"2023-09-22T20:23:11.359642Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Effects of batch size","metadata":{}},{"cell_type":"code","source":"num_genes = 18211\nindices = np.arange(5, 5+num_genes)","metadata":{"execution":{"iopub.status.busy":"2023-09-22T20:14:32.724052Z","iopub.execute_input":"2023-09-22T20:14:32.724553Z","iopub.status.idle":"2023-09-22T20:14:32.731347Z","shell.execute_reply.started":"2023-09-22T20:14:32.724517Z","shell.execute_reply":"2023-09-22T20:14:32.729289Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"batch_sizes = [6, 7, 8, 10, 15, 20, 30, 100, 200, 500, 2000]\nrecords = []\nfor bs in batch_sizes:\n    for run in range(5):\n        shuffle(indices)\n        records.append(\n            [bs, np.mean([mrrmse(filtered_df_train.iloc[:, indices[i*bs: (i+1)*bs]]) for i in range(num_genes // bs)])]\n        )","metadata":{"execution":{"iopub.status.busy":"2023-09-22T20:34:34.420936Z","iopub.execute_input":"2023-09-22T20:34:34.421433Z","iopub.status.idle":"2023-09-22T20:35:59.782957Z","shell.execute_reply.started":"2023-09-22T20:34:34.421398Z","shell.execute_reply":"2023-09-22T20:35:59.781245Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df = pd.DataFrame(records, columns=['batch_size(#genes)', 'MRRMSE'])\nsns.boxenplot(df, x='batch_size(#genes)', y='MRRMSE')","metadata":{"execution":{"iopub.status.busy":"2023-09-22T20:37:21.413831Z","iopub.execute_input":"2023-09-22T20:37:21.414353Z","iopub.status.idle":"2023-09-22T20:37:21.868611Z","shell.execute_reply.started":"2023-09-22T20:37:21.414308Z","shell.execute_reply":"2023-09-22T20:37:21.867387Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Breakdown by compound","metadata":{}},{"cell_type":"code","source":"for compound in compound_list_train:\n    df_by_compound = filtered_df_train[filtered_df_train['sm_name'] == compound]\n    print(f\"{compound}: \\t {mrrmse(df_by_compound)}\")","metadata":{"execution":{"iopub.status.busy":"2023-09-22T19:57:15.174928Z","iopub.execute_input":"2023-09-22T19:57:15.175976Z","iopub.status.idle":"2023-09-22T19:57:15.304865Z","shell.execute_reply.started":"2023-09-22T19:57:15.175940Z","shell.execute_reply":"2023-09-22T19:57:15.303438Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Breakdown by cell type and compound","metadata":{"execution":{"iopub.status.busy":"2023-09-17T18:49:34.864122Z","iopub.execute_input":"2023-09-17T18:49:34.864632Z","iopub.status.idle":"2023-09-17T18:49:34.871974Z","shell.execute_reply.started":"2023-09-17T18:49:34.864590Z","shell.execute_reply":"2023-09-17T18:49:34.870475Z"}}},{"cell_type":"code","source":"for compound in compound_list_train:\n    print(compound)\n\n    for cell_type in ['B cells', 'Myeloid cells']:\n    \n        df_by_compound = filtered_df_train[(filtered_df_train['sm_name'] == compound) * (filtered_df_train['cell_type'] == cell_type)]\n        print(f\"\\t{cell_type}:\\t{mrrmse(df_by_compound)}\")","metadata":{"execution":{"iopub.status.busy":"2023-09-22T19:57:15.307022Z","iopub.execute_input":"2023-09-22T19:57:15.308152Z","iopub.status.idle":"2023-09-22T19:57:15.511720Z","shell.execute_reply.started":"2023-09-22T19:57:15.308098Z","shell.execute_reply":"2023-09-22T19:57:15.510229Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Remove two outlier compounds","metadata":{}},{"cell_type":"code","source":"mrrmse(filtered_df_train[~filtered_df_train['sm_name'].isin(['MLN 2238', 'Dactolisib'])])","metadata":{"execution":{"iopub.status.busy":"2023-09-22T19:57:15.517341Z","iopub.execute_input":"2023-09-22T19:57:15.527780Z","iopub.status.idle":"2023-09-22T19:57:15.568865Z","shell.execute_reply.started":"2023-09-22T19:57:15.527677Z","shell.execute_reply":"2023-09-22T19:57:15.567624Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Some Random EDA","metadata":{}},{"cell_type":"code","source":"import seaborn as sns","metadata":{"execution":{"iopub.status.busy":"2023-09-17T19:53:50.276088Z","iopub.execute_input":"2023-09-17T19:53:50.276926Z","iopub.status.idle":"2023-09-17T19:53:50.284543Z","shell.execute_reply.started":"2023-09-17T19:53:50.276864Z","shell.execute_reply":"2023-09-17T19:53:50.282965Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_four_cells = df_train[df_train[\"cell_type\"].isin([\"NK cells\", \"T cells CD4+\", \"T cells CD8+\", \"T regulatory cells\"])]","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"mrrmse(df_four_cells)","metadata":{"execution":{"iopub.status.busy":"2023-09-17T19:27:32.250875Z","iopub.execute_input":"2023-09-17T19:27:32.251292Z","iopub.status.idle":"2023-09-17T19:27:32.471302Z","shell.execute_reply.started":"2023-09-17T19:27:32.251259Z","shell.execute_reply":"2023-09-17T19:27:32.470165Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for compound in compound_list_train:\n    df_by_compound = df_four_cells[df_four_cells['sm_name'] == compound]\n    print(f\"{compound}: \\t {mrrmse(df_by_compound)}\")","metadata":{"execution":{"iopub.status.busy":"2023-09-17T19:28:06.337023Z","iopub.execute_input":"2023-09-17T19:28:06.337431Z","iopub.status.idle":"2023-09-17T19:28:06.456699Z","shell.execute_reply.started":"2023-09-17T19:28:06.337399Z","shell.execute_reply":"2023-09-17T19:28:06.455411Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"records = []\nfor compound in list(df_train['sm_name'].value_counts().index):\n    df_by_compound = df_four_cells[df_four_cells['sm_name'] == compound]\n    if compound in compound_list_train:\n        records.append([compound, \"train\", mrrmse(df_by_compound)])\n        print(f\"\\t{compound}: \\t {mrrmse(df_by_compound)}\")\n    else:\n        records.append([compound, \"test\", mrrmse(df_by_compound)])\n        print(f\"{compound}: \\t {mrrmse(df_by_compound)}\")","metadata":{"execution":{"iopub.status.busy":"2023-09-17T19:32:46.219237Z","iopub.execute_input":"2023-09-17T19:32:46.219675Z","iopub.status.idle":"2023-09-17T19:32:48.324410Z","shell.execute_reply.started":"2023-09-17T19:32:46.219641Z","shell.execute_reply":"2023-09-17T19:32:48.323018Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"raw","source":"","metadata":{}},{"cell_type":"code","source":"df = pd.DataFrame(records, columns = ['compound', 'train', 'mrrmse'])#.set_index(\"compound\")","metadata":{"execution":{"iopub.status.busy":"2023-09-17T19:37:02.466381Z","iopub.execute_input":"2023-09-17T19:37:02.466873Z","iopub.status.idle":"2023-09-17T19:37:02.474916Z","shell.execute_reply.started":"2023-09-17T19:37:02.466831Z","shell.execute_reply":"2023-09-17T19:37:02.473319Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sns.scatterplot(df, x=\"compound\", y=\"mrrmse\", hue=\"train\")","metadata":{"execution":{"iopub.status.busy":"2023-09-17T19:37:48.830314Z","iopub.execute_input":"2023-09-17T19:37:48.830739Z","iopub.status.idle":"2023-09-17T19:37:50.438859Z","shell.execute_reply.started":"2023-09-17T19:37:48.830704Z","shell.execute_reply":"2023-09-17T19:37:50.437656Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{"execution":{"iopub.status.busy":"2023-09-17T18:52:42.370354Z","iopub.execute_input":"2023-09-17T18:52:42.370808Z","iopub.status.idle":"2023-09-17T18:52:43.778874Z","shell.execute_reply.started":"2023-09-17T18:52:42.370774Z","shell.execute_reply":"2023-09-17T18:52:43.777307Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sns.boxenplot(filtered_df_train.iloc[:10, 5:].T,)","metadata":{"execution":{"iopub.status.busy":"2023-09-17T19:05:23.894727Z","iopub.execute_input":"2023-09-17T19:05:23.895142Z","iopub.status.idle":"2023-09-17T19:05:24.367914Z","shell.execute_reply.started":"2023-09-17T19:05:23.895106Z","shell.execute_reply":"2023-09-17T19:05:24.366388Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sns.boxenplot(filtered_df_train.iloc[10:20, 5:].T,)","metadata":{"execution":{"iopub.status.busy":"2023-09-17T19:05:28.372882Z","iopub.execute_input":"2023-09-17T19:05:28.373313Z","iopub.status.idle":"2023-09-17T19:05:28.825831Z","shell.execute_reply.started":"2023-09-17T19:05:28.373280Z","shell.execute_reply":"2023-09-17T19:05:28.824664Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sns.boxenplot(filtered_df_train.iloc[20:30, 5:].T,)","metadata":{"execution":{"iopub.status.busy":"2023-09-17T19:05:30.837063Z","iopub.execute_input":"2023-09-17T19:05:30.837446Z","iopub.status.idle":"2023-09-17T19:05:31.165076Z","shell.execute_reply.started":"2023-09-17T19:05:30.837415Z","shell.execute_reply":"2023-09-17T19:05:31.163642Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"filtered_df_train.iloc[:, 5:].T.quantile(0.9)","metadata":{"execution":{"iopub.status.busy":"2023-09-17T19:15:45.314763Z","iopub.execute_input":"2023-09-17T19:15:45.315134Z","iopub.status.idle":"2023-09-17T19:15:45.345268Z","shell.execute_reply.started":"2023-09-17T19:15:45.315104Z","shell.execute_reply":"2023-09-17T19:15:45.343922Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}