{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":59094,"databundleVersionId":7010844,"sourceType":"competition"},{"sourceId":7076612,"sourceType":"datasetVersion","datasetId":4070302},{"sourceId":152589679,"sourceType":"kernelVersion"}],"dockerImageVersionId":30587,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"![Human Proteins Atlas](https://user-images.githubusercontent.com/115424463/285138737-dfc99773-a84b-4821-8dd9-1e5c81bfe4d3.jpg)\n## [Immune cells data](https://www.proteinatlas.org/humanproteome/immune+cell)","metadata":{}},{"cell_type":"markdown","source":"__Hello, fellow Kagglers!__\n\n_I was 2 hours late with the publication of the notebook when the last week of the competition began, and no publications were allowed. However, even after the end of the competition, the idea still seems relevant_\n\nI have an idea. The thing is, I joined late, so I won't have time to exhaust it.<br>\nMy concern is that there might be some value for the medical community, both if the assumption is true or strictly false.\n\n`Here is the idea`:\n\n_The task_<br>\nWe have __cell types__ and __compounds__. The researcher wants to know _what will be the genes expression in reaction to a compound, given a certain type of a cell_.\n\n_The proposal_<br>\nLet's consider `cell_type` as a normalization function. We take a vector of genes expression in  response to a particular compound and 'adjust' it based on the specified cell type. After computing this adjustment for each cell type, we normalize each row to a hypothetical \"dummy cell\" condition, eliminating the \"cell type\" information from the data. This data structuring puts us in an optimal position for training and prediction using a single predictor — the compound. Subsequently, the forecasted results of the trained model undergo 'restoration' per cell types, producing the final values.\n\n_Why not standard methods_ <br>\n__Target-transforming__ algorithms (all kinds tree models, KNN, SVM, etc.) do not capture the functional dependence between `X` and `y`, and do not extend beyond the y_values avaliable in the training data. That's not what a researcher is looking for in this case.\n\n__Feature-transforming__ algorithms (such as LinReg, neural networks) do capture the correlation and non-linear dependencies but demand a substantial amount of data for effective training. In our training set, we have 6 cell types, 146 compounds, and 614 rows (to be precise, 514, as 100 rows contain compounds not used in the test set). This results in 87% of compounds being tested on only 4 cell types (see the plot below). When fitting a curve with just 4 data points, the best achievable outcome is the mean (or some quantile, essentially a guess). While attempting to slightly curve the line with a 2-degree polynomial is possible, it largely relies on chance and does not surpass blending (I experimented with it and wasn't successful). Let's not forget about outliers, noise, multicoliniarity among genes, and other distortions. This is why we end up randomly blending zeros with means and quantiles to get a better score. \n\nIn a nutshell, there is not enough CASES (cell_type x compound) for analytical models to derive meaningful insights from the data. The quantity of additional information about cell_types added to the training set doesn't significantly impact the outcomes. For instance, I included gene expression values per cell_type (options: 'nTPM', 'specificity_score') as extra predictors, but it had no discernible effect (although it might have worked with a larger dataset).\n\n`_REMEDY_`<br>\nWe have [immune cells data](https://www.proteinatlas.org/humanproteome/immune+cell) from the Human Protein Atlas. It reads that each cell_type is carecterized by a certain expression of leading genes. Does this imply that the type of a cell influences how its genes respond to any given compound? Essentially, this is the core question of the competition.\n\nSo, here is the thing. I am not a biologist and I have zero domain knowledge. It takes a knowing person  to harmonize the measuring units from the Human Protein Atlas with the data in this competition. Once that's accomplished, we can normalize the data based on cell_type to a hypothetical `dummy_cell` condition and restore it after making predictions. I believe, this can work.    \n\n_Steps taken_ <br>\nAs a normalization factor, I tried several options:\n- Difference between a gene's global mean expression and cell_type-specific mean expression (all genes)\n- The same for the leading genes only (works in the oposite direction)\n- Attempted to enhance the training set by incorporating Atlas data for different models.\n\n_Steps to take_\n- Bring Atlas and competion's data to the same measuring unites and try the `dummy cell normalization` thing\n\n<blockquote style=\"margin-right:auto; margin-left:auto; background-color: #faf0be; padding: 1em; margin:24px;\">\nPlease upvote if you find my ideas worth consideration <br> \n    </blockquote>","metadata":{}},{"cell_type":"markdown","source":"# Imports","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\n\nimport matplotlib.pyplot as plt\nimport seaborn as sns\n\nfrom sklearn.multioutput import MultiOutputRegressor\nfrom sklearn.decomposition import TruncatedSVD\nfrom sklearn.linear_model import LinearRegression\nfrom sklearn.preprocessing import OneHotEncoder\nfrom sklearn.compose import ColumnTransformer\nfrom sklearn.preprocessing import PolynomialFeatures\n\nfrom catboost import CatBoostRegressor, Pool\n\nimport warnings\nwarnings.simplefilter('ignore')\n\nSEED = 34\nnp.random.seed(SEED)","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2023-12-01T10:01:45.094102Z","iopub.execute_input":"2023-12-01T10:01:45.094511Z","iopub.status.idle":"2023-12-01T10:01:45.102404Z","shell.execute_reply.started":"2023-12-01T10:01:45.094483Z","shell.execute_reply":"2023-12-01T10:01:45.101337Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Data","metadata":{}},{"cell_type":"markdown","source":"__Train__:\nAll compounds in T, NK cells\n15 compounds + positive and negative controls in B and myeloid cells\n\n__Test__:\nRandomly selected compounds in B and myeloid cells\n__Note that there is no additional test data beyond the indicated `cell_type / sm_name pairs`__. The input to your model will be a tuple of cell_type and sm_name and the output of your model will be predicted signed -log10(p-values) for all 18211 genes.\n\n__de_train.parquet__ - Aggregated differential expression data in dense array format.\n- genes `A1BG, A1BG-AS1, …, ZZEF1` (numbering __18,211__ in all) - Differential expression value (`-log10(p-value) * sign(LFC)`) for each gene. Here, LFC is the estimated log-fold change in expression between the treatment and control condition after shrinkage as calculated by Limma. Positive LFC means the gene goes up in the treatment condition relative to the control.\n- __`cell_type`__ - The annotated cell type of each cell based on RNA expression.\n-__`sm_name`__ - The primary name for the (parent) compound (in a standardized representation) as chosen by LINCS. This is provided to map the data in this experiment to the LINCS Connectivity Map data.\nsm_lincs_id - The global LINCS ID (parent) compound (in a standardized representation). This is provided to map the data in this experiment to the LINCS Connectivity Map data.\n- SMILES - Simplified molecular-input line-entry system (SMILES) representations of the compounds used in the experiment. This is a 1D representation of molecular structure. These SMILES are provided by Cellarity based on the specific compounds ordered for this experiment.\n","metadata":{}},{"cell_type":"code","source":"train = pd.read_parquet('/kaggle/input/open-problems-single-cell-perturbations/de_train.parquet')\ntest = pd.read_csv('/kaggle/input/open-problems-single-cell-perturbations/id_map.csv')","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:45.104842Z","iopub.execute_input":"2023-12-01T10:01:45.105583Z","iopub.status.idle":"2023-12-01T10:01:46.372853Z","shell.execute_reply.started":"2023-12-01T10:01:45.105531Z","shell.execute_reply":"2023-12-01T10:01:46.371645Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nThe number of compounds tested on a respective number of unique cell types.\n87% are tested on 4 cell types only. \n\"\"\"\ncounts_df = train.groupby('sm_name')['cell_type'].nunique().reset_index(name='cell_type_count')\ncompounds_counts = counts_df['cell_type_count'].value_counts().sort_index()\ncompounds_percentages = compounds_counts / len(counts_df) * 100\nbarplot_df = pd.DataFrame({'Number of Cell Types': compounds_percentages.index,\n                           'Percentage of Compounds': compounds_percentages.values})\nbarplot_df = barplot_df.sort_values(by='Number of Cell Types')\n\nplt.figure(figsize=(6, 5))\n\nax = sns.barplot(x='Number of Cell Types', y='Percentage of Compounds', data=barplot_df, palette='viridis')\n\nfor p in ax.patches:\n    ax.annotate(f'{p.get_height():.2f}%', (p.get_x() + p.get_width() / 2., p.get_height()),\n                ha='center', va='center', xytext=(0, 5), textcoords='offset points', fontsize=11, color='black')\n\nplt.xlabel('Number of Cell Types')\nplt.ylabel('Percentage')\nplt.title('Percentage of Compounds Tested on the Respective Number of Cell Types')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:46.374370Z","iopub.execute_input":"2023-12-01T10:01:46.374715Z","iopub.status.idle":"2023-12-01T10:01:46.614618Z","shell.execute_reply.started":"2023-12-01T10:01:46.374685Z","shell.execute_reply":"2023-12-01T10:01:46.613485Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# train.query('cell_type == \"B cells\" | cell_type == \"Myeloid cells\"').shape\nprint(train.query('cell_type == \"B cells\" | cell_type == \"Myeloid cells\"')['sm_name'].nunique())\nprint(test.query('cell_type == \"B cells\" | cell_type == \"Myeloid cells\"')['sm_name'].nunique())","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:46.617893Z","iopub.execute_input":"2023-12-01T10:01:46.618639Z","iopub.status.idle":"2023-12-01T10:01:48.223662Z","shell.execute_reply.started":"2023-12-01T10:01:46.618595Z","shell.execute_reply":"2023-12-01T10:01:48.222389Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nComparison\n\"\"\"\nprint (f'N compounds in training set {train.sm_name.nunique()}\\n'\n       f'N compounds in test set {test.sm_name.nunique()}')","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:48.225184Z","iopub.execute_input":"2023-12-01T10:01:48.226170Z","iopub.status.idle":"2023-12-01T10:01:48.233516Z","shell.execute_reply.started":"2023-12-01T10:01:48.226107Z","shell.execute_reply":"2023-12-01T10:01:48.232063Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nThese 17 compounds are not present in the test set. \nEssentially, this implies that in the training data, \nwe have 100 fewer rows when excluding these compounds.\n\"\"\"\ndiff = set(train.sm_name.sort_values().unique()) - set(test.sm_name.sort_values().unique())\nprint (len(diff))\nprint (len(train.query('sm_name in @diff')))","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:48.235425Z","iopub.execute_input":"2023-12-01T10:01:48.235790Z","iopub.status.idle":"2023-12-01T10:01:48.857890Z","shell.execute_reply.started":"2023-12-01T10:01:48.235759Z","shell.execute_reply":"2023-12-01T10:01:48.856364Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Dummy Cell approach","metadata":{"execution":{"iopub.status.busy":"2023-11-17T23:21:34.080656Z","iopub.execute_input":"2023-11-17T23:21:34.081570Z","iopub.status.idle":"2023-11-17T23:21:34.092537Z","shell.execute_reply.started":"2023-11-17T23:21:34.081539Z","shell.execute_reply":"2023-11-17T23:21:34.091720Z"}}},{"cell_type":"code","source":"genes = train.columns[5:].to_list()","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:48.859860Z","iopub.execute_input":"2023-12-01T10:01:48.860701Z","iopub.status.idle":"2023-12-01T10:01:48.868436Z","shell.execute_reply.started":"2023-12-01T10:01:48.860656Z","shell.execute_reply":"2023-12-01T10:01:48.867217Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"'''\nThe leading genes identified for each cell type sourced from the Human Protein Atlas.\nCommented are genes that are absent from the competition's data.\n'''\nexpressed_genes = {}\n\nexpressed_genes ['NK cells'] = [\n#       ('COX6A2',15.0),\n#       ('ZMAT4',8.6),\n      ('KIR2DL4', 8.0, 80),\n      ('APOC2',3.2, 29),\n      ('BAALC',3.7, 28),\n      ('SERPINE1',9.5, 25),\n#       ('NFILZ',14.5),\n#       ('EXOC1L',1.9),\n      ('RAMP1', 70.4, 18),\n      ('LGALS9B',13.3, 17),\n#       ('TRIM54',3.6),\n#       ('PRDM16',2.4)\n]\n\nexpressed_genes ['T cells CD4+'] = [\n    ('FHIT',87.9, 4),\n    ('ADAM23',2.1, 16),\n    ('NEFL',11.2, 5)\n]\n\nexpressed_genes['T cells CD8+'] =[\n    ('CFAP97D2',4.0, 13),\n    ('CD248',21.9, 7),\n#     ('MXRA8',6.4),\n#     ('ZCCHC12',1.2),\n    ('NRCAM', 2.1, 6),\n    ('REG4', 24.7, 6),\n    ('SFRP5', 3.6, 4),\n    ('MT3', 5.0, 4)\n]\n\nexpressed_genes ['T regulatory cells'] = [\n#     ('ENSG00000290184', 2.0, 17),\n    ('FANK1', 61.4, 16),\n#     ('SEMA3G', 1.4),\n    ('FOXP3', 52.0, 11),\n    ('HACD1', 55.9, 11),\n#     ('ACTG2', 8.6),\n    ('PMCH', 9.6, 10),\n#     ('BOLL', 1.5),\n    ('TRIM16', 44.9, 9),\n#     ('OGN', 1.7),\n    ('DUSP4', 6.8, 9),\n    ('RTKN2', 16.7, 8)\n]\n\nexpressed_genes ['B cells'] = [\n    ('VPREB3', 807.5, 730),\n#     ('ENSG00000275063', 404.4),\n#     ('ENSG00000277856', 34.5),\n    ('COL19A1', 32.0, 299),\n    ('IGLL5', 619.6, 247),\n#     ('ENSG00000277836', 17.7),\n    ('WNT16', 15.6, 146),\n    ('KLHL14', 10.6, 106),\n    ('MS4A1', 630.0, 96),\n    ('CPNE5', 51.7, 93),\n    ('CD19', 172.5, 90),\n    ('STEAP1B', 8.4, 85)\n]\n\n\nexpressed_genes ['Myeloid cells'] = [\n#     ('XCR1', 2.3),\n    ('MRC1', 4.8, 16),\n    ('CRIP3',44.7, 16),\n    ('CD1E',62.6, 14),\n    ('HTR7',1.4, 12),\n#     ('CD207',1.1),\n    ('CLEC10A',466.5, 12),\n    ('NDRG2',204.5, 11),\n    ('DEPTOR',20.3, 11),\n    ('EHF',1.7, 11),\n    ('PKIB',34.9, 10),\n#     ('C19orf33',22.3)\n]\n\nfor key in expressed_genes.keys():\n    for i in expressed_genes[key]:\n        print(f'{i[0]} is in genes: {i[0] in genes}')","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:48.870032Z","iopub.execute_input":"2023-12-01T10:01:48.870768Z","iopub.status.idle":"2023-12-01T10:01:48.899218Z","shell.execute_reply.started":"2023-12-01T10:01:48.870735Z","shell.execute_reply":"2023-12-01T10:01:48.897117Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"'''\nDictionary to DataFrame\n'''\nexpressed_genes_df = pd.DataFrame.from_dict(expressed_genes, orient = 'index')\nexpressed_genes_df = (expressed_genes_df\n                          .stack()\n                          .to_frame('genes')\n                          .reset_index().drop('level_1', axis=1)\n                          .rename({'level_0':'cell_type'}, axis=1)\n                         )\nexpressed_genes_df [['gene', 'nTPM', 'specificity_score']] = expressed_genes_df['genes'].apply(pd.Series)\nexpressed_genes_df  = expressed_genes_df.drop('genes', axis=1)\n\nexpressed_genes_df.head(3)","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:48.904624Z","iopub.execute_input":"2023-12-01T10:01:48.905102Z","iopub.status.idle":"2023-12-01T10:01:48.939936Z","shell.execute_reply.started":"2023-12-01T10:01:48.905067Z","shell.execute_reply":"2023-12-01T10:01:48.938742Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nLables for expressed genes\n\"\"\"\nall_expressed_genes = expressed_genes_df['gene'].unique()\nall_expressed_genes_with_suffix = [f'{gene}_exp' for gene in all_expressed_genes]","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:48.941632Z","iopub.execute_input":"2023-12-01T10:01:48.942087Z","iopub.status.idle":"2023-12-01T10:01:48.953921Z","shell.execute_reply.started":"2023-12-01T10:01:48.942043Z","shell.execute_reply":"2023-12-01T10:01:48.952787Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nDataFrame for expressed genes in wide format. nTPM only\n\"\"\"\nexpressed = expressed_genes_df.copy().drop('specificity_score', axis=1)\nexpressed = expressed.set_index(['cell_type', 'gene']).unstack('gene').fillna(0)\n\nexpressed.columns = all_expressed_genes_with_suffix\nexpressed.head(2)","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:48.955512Z","iopub.execute_input":"2023-12-01T10:01:48.957811Z","iopub.status.idle":"2023-12-01T10:01:49.015100Z","shell.execute_reply.started":"2023-12-01T10:01:48.957765Z","shell.execute_reply":"2023-12-01T10:01:49.014260Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nDataFrame for expressed genes in wide format. 'Specificity score' only\n\"\"\"\nexpressed_spc = expressed_genes_df.copy().drop('nTPM', axis=1)\nexpressed_spc = expressed_spc.set_index(['cell_type', 'gene']).unstack('gene').fillna(0)\n\nexpressed_spc.columns = all_expressed_genes_with_suffix\nexpressed_spc.head(2)","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:49.017555Z","iopub.execute_input":"2023-12-01T10:01:49.018405Z","iopub.status.idle":"2023-12-01T10:01:49.064553Z","shell.execute_reply.started":"2023-12-01T10:01:49.018360Z","shell.execute_reply":"2023-12-01T10:01:49.063416Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Calculating Normalization Factor \n> as difference between each gene's global average and cell-spesific mean ","metadata":{}},{"cell_type":"code","source":"cell_specific_mean = train.groupby('cell_type')[genes].mean()\nglobal_avg = train[genes].mean().to_frame().T\nnorm_factor = (-1) * cell_specific_mean.sub(global_avg.iloc[0], axis=1)\nnorm_factor.head(3)","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:49.066338Z","iopub.execute_input":"2023-12-01T10:01:49.067020Z","iopub.status.idle":"2023-12-01T10:01:49.392890Z","shell.execute_reply.started":"2023-12-01T10:01:49.066981Z","shell.execute_reply":"2023-12-01T10:01:49.392068Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nNormalization factor for all genes (long format)\n\"\"\"\nnorm_factor_long = (norm_factor\n                    .stack()\n                    .to_frame()\n                    .reset_index()\n                    .rename({'level_1':'gene', 0:'norm_factor'}, axis=1)\n                   )\nnorm_factor_long.head(4)","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:49.393965Z","iopub.execute_input":"2023-12-01T10:01:49.394472Z","iopub.status.idle":"2023-12-01T10:01:49.423176Z","shell.execute_reply.started":"2023-12-01T10:01:49.394441Z","shell.execute_reply":"2023-12-01T10:01:49.422045Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"'''\nNormalization factor for expressed genes (wide format)\n'''\nleft_ = expressed_genes_df[['cell_type','gene']]\n\nexpressed_genes_norm = (pd.merge(left_, norm_factor_long, how='left')\n                      .set_index(['cell_type','gene'])\n                      .reset_index(['cell_type','gene'])\n                      .pivot(index='cell_type', columns = 'gene')\n                      .fillna(0)\n                     )\nexpressed_genes_norm.columns = all_expressed_genes\n\n\nnorm_factor_broad  = pd.DataFrame(0, index = test.set_index(['id','cell_type']).index, columns = genes).reset_index('id')\nnorm_factor_broad.update(expressed_genes_norm)\nnorm_factor_broad = norm_factor_broad.reset_index().set_index(['id', 'cell_type'])\nnorm_factor_broad[all_expressed_genes].head(3)","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:49.424633Z","iopub.execute_input":"2023-12-01T10:01:49.425077Z","iopub.status.idle":"2023-12-01T10:01:49.644201Z","shell.execute_reply.started":"2023-12-01T10:01:49.425047Z","shell.execute_reply":"2023-12-01T10:01:49.643048Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Visualization","metadata":{}},{"cell_type":"code","source":"left_ = expressed_genes_df[['cell_type','gene']]\n\nright_1 = global_avg.T.reset_index().rename({'index':'gene', 0: \"gl_avg\"}, axis =1)\nviz_data = pd.merge(left_, right_1, on ='gene', how = 'left')\n\nright_2 = cell_specific_mean.stack().to_frame('cell_specific_avg').reset_index().rename({'level_1':'gene'}, axis = 1)\nviz_data = pd.merge(viz_data, right_2, on =['cell_type','gene'], how = 'left')\n\ndisplay (viz_data.shape)\nviz_data.head(2)","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:49.646062Z","iopub.execute_input":"2023-12-01T10:01:49.647268Z","iopub.status.idle":"2023-12-01T10:01:49.753350Z","shell.execute_reply.started":"2023-12-01T10:01:49.647227Z","shell.execute_reply":"2023-12-01T10:01:49.752551Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for cell in train['cell_type'].unique():\n    plt.figure(figsize=(8, 4))\n    plt.xticks(rotation=25)\n    plt.title(f'{cell} and its expressed genes against global average')\n\n    data_ = viz_data.query('cell_type == @cell')\n\n    ax = sns.barplot(data_, x='gene', y='gl_avg', color='green', label='global_avg', dodge=True)\n    sns.barplot(data_, x='gene', y='cell_specific_avg',color='yellow', label='cell_specific_avg', alpha=0.3, dodge = True, ax=ax)\n    ax.legend(loc='upper right', fontsize='small')\n    for p in ax.patches:\n            ax.annotate(f'{p.get_height():.2f}', (p.get_x() + p.get_width() / 2., p.get_height()),\n                        ha='center', va='center', xytext=(0, 10), textcoords='offset points', fontsize=8, color='black')\n    \n    ax.set_ylim(ax.get_ylim()[0], ax.get_ylim()[1] * 1.3)\n    ax.set_ylabel(None)\n    ax.set_xlabel(None)\n    \n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:49.754441Z","iopub.execute_input":"2023-12-01T10:01:49.755019Z","iopub.status.idle":"2023-12-01T10:01:51.885230Z","shell.execute_reply.started":"2023-12-01T10:01:49.754990Z","shell.execute_reply":"2023-12-01T10:01:51.884122Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Baseline\n> Publicly avalible blend to serve as a baseline to test the hypothesis <br>\n> Score: 0.567","metadata":{}},{"cell_type":"code","source":"\"\"\" \nPublicly available blends\n\"\"\"\ndf = pd.read_csv('/kaggle/input/op-scp-competition/submission_0531.csv')\ndf_blend = test[['cell_type']].join(df).set_index(['id','cell_type'])","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:51.886700Z","iopub.execute_input":"2023-12-01T10:01:51.887072Z","iopub.status.idle":"2023-12-01T10:01:58.518862Z","shell.execute_reply.started":"2023-12-01T10:01:51.887042Z","shell.execute_reply":"2023-12-01T10:01:58.517599Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Scheme_1: Normalization factor to all genes\n> Adding  - Score: 1.01 <br>\n> Substracting - Score: 0.937","metadata":{}},{"cell_type":"code","source":"cell_types = norm_factor.index.to_list()\ndfs = {}\n\nfor cell in cell_types:\n    dfs[cell] = df_blend.query('cell_type == @cell')\n    dfs[cell] = dfs[cell].add(norm_factor.loc[cell], axis=1)\n    dfs[cell] = dfs[cell].sub(norm_factor.loc[cell], axis=1)\n\ndf_final = pd.concat(dfs.values(), axis=0).reset_index(['cell_type'], drop=True)\n\n# df_final.to_csv('submission.csv')\ndf_final.head(4)","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:01:58.520362Z","iopub.execute_input":"2023-12-01T10:01:58.520726Z","iopub.status.idle":"2023-12-01T10:02:03.114718Z","shell.execute_reply.started":"2023-12-01T10:01:58.520694Z","shell.execute_reply":"2023-12-01T10:02:03.113503Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Scheme_2: Normalization factor to expressed genes only\n> Score: 0.568","metadata":{}},{"cell_type":"code","source":"df_final = df_blend + norm_factor_broad\n\nsubmission = df_final.reset_index(['cell_type'], drop=True)\n# submission.to_csv('submission.csv')","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:02:03.116552Z","iopub.execute_input":"2023-12-01T10:02:03.117036Z","iopub.status.idle":"2023-12-01T10:02:16.035203Z","shell.execute_reply.started":"2023-12-01T10:02:03.116991Z","shell.execute_reply":"2023-12-01T10:02:16.033937Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Scheme_3: Comparizon with 'all zeros' baseline\n> Baseline: 0.666 <br>\n> Score: 0.667","metadata":{}},{"cell_type":"code","source":"baseline = pd.read_csv('/kaggle/input/open-problems-single-cell-perturbations/sample_submission.csv')\n\nbaseline = test[['cell_type']].join(baseline).set_index(['id','cell_type'])\n\ndf_final = baseline + norm_factor_broad\n\nsubmission = df_final.reset_index(['cell_type'], drop=True)\n# submission.to_csv('submission.csv')","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:02:16.036614Z","iopub.execute_input":"2023-12-01T10:02:16.037050Z","iopub.status.idle":"2023-12-01T10:02:21.082714Z","shell.execute_reply.started":"2023-12-01T10:02:16.037008Z","shell.execute_reply":"2023-12-01T10:02:21.081634Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Scheme_4: Teaching with expressed genes per cell_type (nTPM)\n> Baseline (mean): 0.663 <br>\n> Score: 0.92","metadata":{}},{"cell_type":"code","source":"train_extended = pd.merge(train, expressed, on='cell_type', how='left')\ntest_extended = pd.merge(test, expressed, on='cell_type', how='left')\n\n\nmodel = MultiOutputRegressor(\n    CatBoostRegressor(\n        cat_features = ['cell_type', 'sm_name'],\n        verbose = 0,\n        iterations = 200,\n        depth = 4,\n        learning_rate=0.01,\n        loss_function = 'RMSE'\n    )\n)\nX_train = train_extended [['cell_type', 'sm_name'] + all_expressed_genes_with_suffix]\ny_train = train_extended .loc[:,genes].values\n\nreducer = TruncatedSVD (n_components = 100, n_iter = 7, random_state = 34)\ny_train_r = reducer.fit_transform(y_train)\n\nmodel.fit(X_train, y_train_r)\n\nX_test = test_extended[['cell_type', 'sm_name']+ all_expressed_genes_with_suffix]\ny_test_r = model.predict(X_test)\ny_test = reducer.inverse_transform(y_test_r)\n\ndf_test = pd.DataFrame(y_test, columns = genes)\ndf_test.index.name = 'id'\n# df_test.to_csv('submission.csv')","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:02:21.084384Z","iopub.execute_input":"2023-12-01T10:02:21.085586Z","iopub.status.idle":"2023-12-01T10:02:46.892643Z","shell.execute_reply.started":"2023-12-01T10:02:21.085538Z","shell.execute_reply":"2023-12-01T10:02:46.890504Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Scheme_5: Teaching with expressed genes per cell_type (specificity_score)\n> Baseline (mean): 0.663 <br>\n> Score: 0.832","metadata":{}},{"cell_type":"code","source":"train_extended = pd.merge(train, expressed_spc, on='cell_type', how='left')\ntest_extended = pd.merge(test, expressed_spc, on='cell_type', how='left')\n\nX_train = train_extended [['cell_type', 'sm_name'] + all_expressed_genes_with_suffix]\ny_train = train_extended .loc[:,genes].values\n\nreducer = TruncatedSVD (n_components = 100, n_iter = 7, random_state = 34)\ny_train_r = reducer.fit_transform(y_train)\n\nmodel.fit(X_train, y_train_r)\n\nX_test = test_extended[['cell_type', 'sm_name']+ all_expressed_genes_with_suffix]\ny_test_r = model.predict(X_test)\ny_test = reducer.inverse_transform(y_test_r)\n\ndf_test = pd.DataFrame(y_test, columns = genes)\ndf_test.index.name = 'id'\n# df_test.to_csv('submission.csv')","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:02:46.895590Z","iopub.execute_input":"2023-12-01T10:02:46.896788Z","iopub.status.idle":"2023-12-01T10:03:13.320841Z","shell.execute_reply.started":"2023-12-01T10:02:46.896726Z","shell.execute_reply":"2023-12-01T10:03:13.319204Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Scheme_6: Teaching with expressed genes per cell_type (specificity_score) with LinReg ()\n> Baseline (mean): 0.663 <br>\n> Score: 0.779","metadata":{}},{"cell_type":"code","source":"train_extended = pd.merge(train, expressed, on='cell_type', how='left')\ntest_extended = pd.merge(test, expressed, on='cell_type', how='left')\n\nX_train = train_extended[['cell_type', 'sm_name'] + all_expressed_genes_with_suffix]\ny_train = train_extended.loc[:, genes].values\n\nX_test = test_extended[['cell_type', 'sm_name'] + all_expressed_genes_with_suffix]\n\ncategorical_features = ['cell_type', 'sm_name']\ncategorical_transformer = OneHotEncoder(handle_unknown='ignore')\n\npreprocessor = ColumnTransformer(\n    transformers=[\n        ('cat', categorical_transformer, categorical_features)\n    ])\n\nX_train_transformed = preprocessor.fit_transform(X_train)\nreducer = TruncatedSVD(n_components=100, n_iter=7, random_state=34)\ny_train_r = reducer.fit_transform(y_train)\n\nlr = LinearRegression()\n\nlr.fit(X_train_transformed, y_train_r)\nX_test_transformed = preprocessor.transform(X_test)\ny_test_r = lr.predict(X_test_transformed)\ny_test = reducer.inverse_transform(y_test_r)\n\ndf_test = pd.DataFrame(y_test, columns=genes)\ndf_test.index.name = 'id'\n# df_test.to_csv('submission.csv')","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:03:13.327758Z","iopub.execute_input":"2023-12-01T10:03:13.328901Z","iopub.status.idle":"2023-12-01T10:03:17.314820Z","shell.execute_reply.started":"2023-12-01T10:03:13.328854Z","shell.execute_reply":"2023-12-01T10:03:17.313246Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Scheme_7: Teaching with expressed genes per cell_type (specificity_score) with LinReg() and 2-degree polinomial\n\n> Baseline (mean): 0.663 <br>\n> Score: 0.784","metadata":{}},{"cell_type":"code","source":"train_extended = pd.merge(train, expressed, on='cell_type', how='left')\ntest_extended = pd.merge(test, expressed, on='cell_type', how='left')\n\nX_train = train_extended[['cell_type', 'sm_name'] + all_expressed_genes_with_suffix]\ny_train = train_extended.loc[:, genes].values\n\nX_test = test_extended[['cell_type', 'sm_name'] + all_expressed_genes_with_suffix]\n\ncategorical_features = ['cell_type', 'sm_name']\ncategorical_transformer = OneHotEncoder(handle_unknown='ignore')\n\npreprocessor = ColumnTransformer(\n    transformers=[\n        ('cat', categorical_transformer, categorical_features)\n    ])\n\nX_train_transformed = preprocessor.fit_transform(X_train)\n\npoly = PolynomialFeatures(degree=2)\nX_train_transformed_poly = poly.fit_transform(X_train_transformed)\n\nreducer = TruncatedSVD(n_components=100, n_iter=7, random_state=34)\ny_train_r = reducer.fit_transform(y_train)\n\nlr = LinearRegression()\n\nlr.fit(X_train_transformed_poly, y_train_r)\n\nX_test_transformed = preprocessor.transform(X_test)\n\n# Adding polynomial features (quadratic)\nX_test_transformed_poly = poly.transform(X_test_transformed)\n\ny_test_r = lr.predict(X_test_transformed_poly)\ny_test = reducer.inverse_transform(y_test_r)\n\ndf_test = pd.DataFrame(y_test, columns=genes)\ndf_test.index.name = 'id'\n# df_test.to_csv('submission.csv')","metadata":{"execution":{"iopub.status.busy":"2023-12-01T10:03:17.317510Z","iopub.execute_input":"2023-12-01T10:03:17.318743Z","iopub.status.idle":"2023-12-01T10:03:21.861357Z","shell.execute_reply.started":"2023-12-01T10:03:17.318678Z","shell.execute_reply":"2023-12-01T10:03:21.859789Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<blockquote style=\"margin-right:auto; margin-left:auto; background-color: #faf0be; padding: 1em; margin:24px;\">\nPlease upvote if you find my ideas worthy of consideration and I would be really glad to learn your thoughts about this idea of mine in the comments.🤝<br> \n    </blockquote>","metadata":{"execution":{"iopub.status.busy":"2023-11-24T21:22:45.625854Z","iopub.execute_input":"2023-11-24T21:22:45.626300Z","iopub.status.idle":"2023-11-24T21:22:45.660371Z","shell.execute_reply.started":"2023-11-24T21:22:45.626267Z","shell.execute_reply":"2023-11-24T21:22:45.658765Z"}}}]}