{"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":"# OP2 - Genes ontology analysis","metadata":{}},{"cell_type":"markdown","source":"In this notebook, I perform gene ontology analysis on the differential expression data from the 'Open Problems - Single-Cell Perturbations' competition, which aims to predict how small molecule perturbations affect gene expression in different cell types.\nIn my analysis, I filter the differential expression p-values at 0.05 significance and identify the enriched genes for each cell type and compound treatment. I then analyze the gene lists using Enrichr to identify enriched gene ontology terms and pathways from the KEGG database.\nI combine the results into a DataFrame containing the enriched terms, adjusted p-values, and associated genes for each cell type and compound pair.\n","metadata":{}},{"cell_type":"code","source":"!pip install --quiet gseapy","metadata":{"execution":{"iopub.status.busy":"2023-10-21T19:03:02.974201Z","iopub.execute_input":"2023-10-21T19:03:02.974578Z","iopub.status.idle":"2023-10-21T19:03:17.988856Z","shell.execute_reply.started":"2023-10-21T19:03:02.974548Z","shell.execute_reply":"2023-10-21T19:03:17.987416Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Load libraries\nimport datetime\nimport pandas as pd\nimport gseapy as gp","metadata":{"execution":{"iopub.status.busy":"2023-10-21T19:03:26.682562Z","iopub.execute_input":"2023-10-21T19:03:26.683006Z","iopub.status.idle":"2023-10-21T19:03:27.308951Z","shell.execute_reply.started":"2023-10-21T19:03:26.682966Z","shell.execute_reply":"2023-10-21T19:03:27.307940Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Load data\nde_train =   pd.read_parquet(\"/kaggle/input/open-problems-single-cell-perturbations/de_train.parquet\")\nid_map = pd.read_csv(\"/kaggle/input/open-problems-single-cell-perturbations/id_map.csv\")","metadata":{"execution":{"iopub.status.busy":"2023-10-21T19:04:00.509318Z","iopub.execute_input":"2023-10-21T19:04:00.509711Z","iopub.status.idle":"2023-10-21T19:04:03.110062Z","shell.execute_reply.started":"2023-10-21T19:04:00.509681Z","shell.execute_reply":"2023-10-21T19:04:03.109015Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# pivot to long\nid_vars = ['cell_type', 'sm_name', 'sm_lincs_id', 'SMILES', 'control']\ndf_long = pd.melt(de_train, id_vars=id_vars, var_name='Gene', value_name='p_value')\ndf_long.head()","metadata":{"execution":{"iopub.status.busy":"2023-10-21T19:05:35.050636Z","iopub.execute_input":"2023-10-21T19:05:35.051026Z","iopub.status.idle":"2023-10-21T19:05:40.533866Z","shell.execute_reply.started":"2023-10-21T19:05:35.050996Z","shell.execute_reply":"2023-10-21T19:05:40.532722Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"grouped = df_long.groupby(['cell_type', 'SMILES'])['Gene'].apply(list).reset_index()\ngrouped.head()","metadata":{"execution":{"iopub.status.busy":"2023-10-21T19:06:47.450245Z","iopub.execute_input":"2023-10-21T19:06:47.450792Z","iopub.status.idle":"2023-10-21T19:06:51.203954Z","shell.execute_reply.started":"2023-10-21T19:06:47.450748Z","shell.execute_reply":"2023-10-21T19:06:51.202585Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"grouped_genes = grouped.groupby(['cell_type', 'SMILES'])['Gene'].sum()","metadata":{"execution":{"iopub.status.busy":"2023-10-21T19:13:41.143488Z","iopub.execute_input":"2023-10-21T19:13:41.143900Z","iopub.status.idle":"2023-10-21T19:13:41.154273Z","shell.execute_reply.started":"2023-10-21T19:13:41.143870Z","shell.execute_reply":"2023-10-21T19:13:41.153375Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"grouped_genes.head()","metadata":{"execution":{"iopub.status.busy":"2023-10-21T19:13:48.781230Z","iopub.execute_input":"2023-10-21T19:13:48.781615Z","iopub.status.idle":"2023-10-21T19:13:48.793678Z","shell.execute_reply.started":"2023-10-21T19:13:48.781588Z","shell.execute_reply":"2023-10-21T19:13:48.792393Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"results = []\n\nfor (cell_type, smiles), genes in grouped_genes.items():\n    enrichment_result = gp.enrichr(gene_list=genes,  \n                                   gene_sets='KEGG_2019_Human',  \n                                   no_plot=True,\n                                   cutoff=0.01)\n    \n    # Add cell_type and SMILES columns to the result DataFrame\n    enrichment_result.res2d['cell_type'] = cell_type\n    enrichment_result.res2d['SMILES'] = smiles\n\n    results.append(enrichment_result.res2d)\n\n# Concatenate results\ncombined_results = pd.concat(results, axis=0)\n\n# Reorder columns to have 'cell_type' and 'SMILES' as the first columns\ncols = ['cell_type', 'SMILES'] + [col for col in combined_results if col not in ['cell_type', 'SMILES']]\ncombined_results = combined_results[cols]","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"combined_results","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"combined_results.to_csv(\"/kaggle/working/genes_ontology.csv\")","metadata":{},"execution_count":null,"outputs":[]}]}