{"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-16T09:49:00.795679Z","iopub.execute_input":"2023-09-16T09:49:00.796105Z","iopub.status.idle":"2023-09-16T09:49:01.438657Z","shell.execute_reply.started":"2023-09-16T09:49:00.796056Z","shell.execute_reply":"2023-09-16T09:49:01.437325Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# # # **IMPORT LIBRARY**","metadata":{}},{"cell_type":"code","source":"%pip install rdkit","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:49:01.441266Z","iopub.execute_input":"2023-09-16T09:49:01.441986Z","iopub.status.idle":"2023-09-16T09:49:15.247718Z","shell.execute_reply.started":"2023-09-16T09:49:01.441923Z","shell.execute_reply":"2023-09-16T09:49:15.245855Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# # # **Cell Perturbations: Machine Learning Analysis**","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nfrom rdkit import Chem\nfrom rdkit.Chem import Draw\nfrom rdkit.Chem import PandasTools\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.linear_model import LinearRegression\nfrom sklearn.metrics import mean_squared_error\nfrom sklearn.preprocessing import OneHotEncoder\nfrom sklearn.preprocessing import StandardScaler\n\n\nimport warnings\nwarnings.filterwarnings(\"ignore\")\n\nimport math\n\nrc = {\n    \"axes.facecolor\": \"#E6FFE6\",\n    \"figure.facecolor\": \"#FFFFFF\",\n    \"axes.edgecolor\": \"#000000\",\n    \"grid.color\": \"#EBEBE7\",\n    \"font.family\": \"serif\",\n    \"axes.labelcolor\": \"#000000\",\n    \"xtick.color\": \"#000000\",\n    \"ytick.color\": \"#000000\",\n    \"grid.alpha\": 0.4\n}\n\nsns.set(rc=rc)\n\nfrom colorama import Style, Fore\nred = Style.BRIGHT + Fore.RED\nblu = Style.BRIGHT + Fore.BLUE\nmgt = Style.BRIGHT + Fore.MAGENTA\ngld = Style.BRIGHT + Fore.YELLOW\nres = Style.RESET_ALL","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:49:15.249876Z","iopub.execute_input":"2023-09-16T09:49:15.250329Z","iopub.status.idle":"2023-09-16T09:49:18.112022Z","shell.execute_reply.started":"2023-09-16T09:49:15.250289Z","shell.execute_reply":"2023-09-16T09:49:18.110822Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"d_train = pd.read_parquet('/kaggle/input/open-problems-single-cell-perturbations/de_train.parquet')\nd_train.head()","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:49:18.113605Z","iopub.execute_input":"2023-09-16T09:49:18.113990Z","iopub.status.idle":"2023-09-16T09:49:21.027253Z","shell.execute_reply.started":"2023-09-16T09:49:18.113946Z","shell.execute_reply":"2023-09-16T09:49:21.026163Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"da_train = pd.read_parquet('/kaggle/input/open-problems-single-cell-perturbations/adata_train.parquet')\nda_train.head()","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:49:21.029979Z","iopub.execute_input":"2023-09-16T09:49:21.030613Z","iopub.status.idle":"2023-09-16T09:50:35.650349Z","shell.execute_reply.started":"2023-09-16T09:49:21.030571Z","shell.execute_reply":"2023-09-16T09:50:35.648905Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"id_map = pd.read_csv('/kaggle/input/open-problems-single-cell-perturbations/id_map.csv')\nid_map.head()","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:50:35.651802Z","iopub.execute_input":"2023-09-16T09:50:35.652326Z","iopub.status.idle":"2023-09-16T09:50:35.688695Z","shell.execute_reply.started":"2023-09-16T09:50:35.652280Z","shell.execute_reply":"2023-09-16T09:50:35.687358Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"submission_df = pd.read_csv('/kaggle/input/open-problems-single-cell-perturbations/sample_submission.csv')\nsubmission_df.head()","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:50:35.690762Z","iopub.execute_input":"2023-09-16T09:50:35.691337Z","iopub.status.idle":"2023-09-16T09:50:39.733547Z","shell.execute_reply.started":"2023-09-16T09:50:35.691286Z","shell.execute_reply":"2023-09-16T09:50:39.732554Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Check for missing values in id_map\nid_map_missing = id_map.isnull().sum()\nprint(\"Missing values in id_map:\")\nprint(id_map_missing)","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:50:39.735156Z","iopub.execute_input":"2023-09-16T09:50:39.735579Z","iopub.status.idle":"2023-09-16T09:50:39.744786Z","shell.execute_reply.started":"2023-09-16T09:50:39.735540Z","shell.execute_reply":"2023-09-16T09:50:39.742522Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Check for missing values in data\nd_train_missing = d_train.isnull().sum()\nprint(\"\\nMissing values in d_train:\")\nprint(d_train_missing)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Drop duplicate rows based on all columns\nd_train.drop_duplicates(inplace=True)","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:50:39.746663Z","iopub.execute_input":"2023-09-16T09:50:39.747081Z","iopub.status.idle":"2023-09-16T09:50:44.012631Z","shell.execute_reply.started":"2023-09-16T09:50:39.747011Z","shell.execute_reply":"2023-09-16T09:50:44.011226Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"d_train.describe().style.background_gradient(cmap='tab20c')","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:50:44.014389Z","iopub.execute_input":"2023-09-16T09:50:44.014936Z","iopub.status.idle":"2023-09-16T09:51:43.488967Z","shell.execute_reply.started":"2023-09-16T09:50:44.014884Z","shell.execute_reply":"2023-09-16T09:51:43.485739Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# # # Exploratory Data Analysis (EDA)","metadata":{}},{"cell_type":"code","source":"# Count the occurrences of each cell type\ncell_type_counts = d_train['cell_type'].value_counts()\n\n# Create a pie chart\nplt.figure(figsize=(10, 6))\nplt.pie(cell_type_counts, labels=cell_type_counts.index, autopct='%1.1f%%', startangle=140, colors=['#ff9999','#66b3ff','#99ff99','#ffcc99','#c2c2f0'])\nplt.title('Distribution of Cell Types', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\nplt.axis('equal')  # Equal aspect ratio ensures the pie chart is circular.\nplt.savefig('Distribution of Cell Types.png')\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:51:43.491638Z","iopub.execute_input":"2023-09-16T09:51:43.492504Z","iopub.status.idle":"2023-09-16T09:51:43.863021Z","shell.execute_reply.started":"2023-09-16T09:51:43.492417Z","shell.execute_reply":"2023-09-16T09:51:43.861547Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(10, 6))\nsns.histplot(data=d_train, x='A1BG', hue='cell_type', kde=True)\nplt.title('Distribution of A1BG Gene Expression', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\nplt.savefig('Distribution of A1BG Gene Expression.png')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:51:43.865228Z","iopub.execute_input":"2023-09-16T09:51:43.867311Z","iopub.status.idle":"2023-09-16T09:51:46.813256Z","shell.execute_reply.started":"2023-09-16T09:51:43.867237Z","shell.execute_reply":"2023-09-16T09:51:46.811706Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Function to convert SMILES to RDKit molecule\ndef smiles_to_mol(smiles):\n    return Chem.MolFromSmiles(smiles)\n\n# Function to draw molecule from RDKit molecule object\ndef draw_molecule(mol, title, save_path):\n    img = Draw.MolToImage(mol)\n    img = img.resize((400, 400))  # Resize to a reasonable size\n    img.save(save_path)  # Save the image\n    plt.figure(figsize=(4, 4))\n    plt.imshow(img)\n    plt.title(title)\n    plt.axis('off')\n    plt.show()\n\n\n# Randomly select 10 different images\nsampled_df = d_train.sample(n=10, random_state=42)\n\n# Example usage of the functions\nfor i, row in sampled_df.iterrows():\n    mol = smiles_to_mol(row['SMILES'])\n    if mol:\n        image_path = f'{row[\"sm_name\"]}.png'  # Define the image path\n        draw_molecule(mol, row['sm_name'], image_path)","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:51:46.815024Z","iopub.execute_input":"2023-09-16T09:51:46.815424Z","iopub.status.idle":"2023-09-16T09:51:49.693639Z","shell.execute_reply.started":"2023-09-16T09:51:46.815384Z","shell.execute_reply":"2023-09-16T09:51:49.692202Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Create a heatmap\nplt.figure(figsize=(8, 7))\nsns.heatmap(d_train.drop(['cell_type', 'sm_name', 'sm_lincs_id', 'SMILES', 'control'], axis=1))\nplt.title('Gene Expression Heatmap', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\nplt.savefig('Gene Expression Heatmap.png')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:51:49.699040Z","iopub.execute_input":"2023-09-16T09:51:49.699443Z","iopub.status.idle":"2023-09-16T09:52:15.259617Z","shell.execute_reply.started":"2023-09-16T09:51:49.699409Z","shell.execute_reply":"2023-09-16T09:52:15.257548Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Create a scatter plot\nplt.figure(figsize=(10, 6))\nsns.scatterplot(x='A1BG', y='A2M', data=d_train, hue='cell_type')\nplt.title('A1BG vs A2M Gene Expression', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\nplt.savefig('A1BG vs A2M Gene Expression.png')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:52:15.262110Z","iopub.execute_input":"2023-09-16T09:52:15.262633Z","iopub.status.idle":"2023-09-16T09:52:16.229331Z","shell.execute_reply.started":"2023-09-16T09:52:15.262582Z","shell.execute_reply":"2023-09-16T09:52:16.227974Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Create a box plot\nplt.figure(figsize=(10, 6))\nsns.boxplot(x='cell_type', y='A1BG', data=d_train)\nplt.title('A1BG Gene Expression Across Cell Types', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\nplt.xticks(rotation=45)\nplt.savefig('A1BG Gene Expression Across Cell Types.png')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:52:16.231272Z","iopub.execute_input":"2023-09-16T09:52:16.232074Z","iopub.status.idle":"2023-09-16T09:52:16.724796Z","shell.execute_reply.started":"2023-09-16T09:52:16.232026Z","shell.execute_reply":"2023-09-16T09:52:16.723327Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Create a violin plot\nplt.figure(figsize=(10, 6))\nsns.violinplot(x='cell_type', y='A2M', data=d_train)\nplt.title('A2M Gene Expression Across Cell Types', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\nplt.xticks(rotation=45)\nplt.savefig('A2M Gene Expression Across Cell Types.png')\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:52:16.726565Z","iopub.execute_input":"2023-09-16T09:52:16.727443Z","iopub.status.idle":"2023-09-16T09:52:17.440280Z","shell.execute_reply.started":"2023-09-16T09:52:16.727391Z","shell.execute_reply":"2023-09-16T09:52:17.438961Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sns.pairplot(d_train[['A1BG', 'A2M', 'A2M-AS1', 'A2MP1']])\nplt.suptitle('Pair Plot of Gene Expressions', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:52:17.441657Z","iopub.execute_input":"2023-09-16T09:52:17.442043Z","iopub.status.idle":"2023-09-16T09:52:24.541198Z","shell.execute_reply.started":"2023-09-16T09:52:17.441978Z","shell.execute_reply":"2023-09-16T09:52:24.539891Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(6, 4))\nsns.heatmap(d_train[['A1BG', 'A2M', 'A2M-AS1', 'A2MP1']].corr(), annot=True, cmap='coolwarm')\nplt.title('Correlation Heatmap', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\nplt.savefig('Correlation Heatmap.png')\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:52:24.543574Z","iopub.execute_input":"2023-09-16T09:52:24.544088Z","iopub.status.idle":"2023-09-16T09:52:25.093242Z","shell.execute_reply.started":"2023-09-16T09:52:24.544036Z","shell.execute_reply":"2023-09-16T09:52:25.091908Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(8, 6))\nsns.clustermap(d_train.drop(['cell_type', 'sm_name', 'sm_lincs_id', 'SMILES', 'control'], axis=1), cmap='coolwarm')\nplt.title('Gene Expression Clustermap', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\nplt.savefig('Gene Expression Clustermap.png')\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:52:25.095081Z","iopub.execute_input":"2023-09-16T09:52:25.095644Z","iopub.status.idle":"2023-09-16T09:54:48.778766Z","shell.execute_reply.started":"2023-09-16T09:52:25.095593Z","shell.execute_reply":"2023-09-16T09:54:48.777437Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from rdkit import Chem\n\n# Define the substructure pattern using SMARTS notation\npatt = Chem.MolFromSmarts('ClccccF')\n\n# The `d_train` is your pandas DataFrame containing SMILES strings in a column named 'SMILES'\nsubms = []\n\nfor smi in d_train['SMILES']:\n    mol = Chem.MolFromSmiles(smi)\n    if mol:\n        hit_ats = list(mol.GetSubstructMatch(patt))\n        if hit_ats:\n            subms.append(mol)\n\n# Print the number of molecules with the substructure\nprint(len(subms))\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:54:48.780368Z","iopub.execute_input":"2023-09-16T09:54:48.780854Z","iopub.status.idle":"2023-09-16T09:54:48.955255Z","shell.execute_reply.started":"2023-09-16T09:54:48.780814Z","shell.execute_reply":"2023-09-16T09:54:48.953773Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# **Chemical Reaction Computations**\n\nUse RDKit to perform a simple reaction computation using SMILES ","metadata":{}},{"cell_type":"code","source":"from rdkit import Chem\n\n# Define a function to perform a chemical reaction\ndef perform_reaction(reactant_smiles, product_smiles, reaction_type):\n    # Convert SMILES strings to RDKit molecules\n    reactant = Chem.MolFromSmiles(reactant_smiles)\n    product = Chem.MolFromSmiles(product_smiles)\n    \n    # This is just an example reaction (dummy reaction for illustration purposes)\n    reaction = Chem.CombineMols(reactant, product)\n    \n    return Chem.MolToSmiles(reaction)\n\n# Example usage\nreactant_smiles = \"CC(=O)OC1=CC=CC=C1C(=O)O\"\nproduct_smiles = \"CC(=O)OC1=CC=CC=C1C(=O)OC\"\n\nresult_smiles = perform_reaction(reactant_smiles, product_smiles, 'example_reaction')\nprint(f\"Resulting SMILES string: {result_smiles}\")","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:54:48.957269Z","iopub.execute_input":"2023-09-16T09:54:48.957826Z","iopub.status.idle":"2023-09-16T09:54:48.967525Z","shell.execute_reply.started":"2023-09-16T09:54:48.957770Z","shell.execute_reply":"2023-09-16T09:54:48.966045Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# **Chemical Property Calculation**\n\nCompute various chemical properties like molecular weight, logP, etc., for each compound","metadata":{}},{"cell_type":"code","source":"from rdkit import Chem\nfrom rdkit.Chem import Descriptors\nimport pandas as pd\n\n# DataFrame 'd_train with a column 'SMILES'\nd_train['Molecule'] = d_train['SMILES'].apply(Chem.MolFromSmiles)\n\n# Check if the 'Molecule' column is correctly created\nprint(d_train.head())\n\n# Calculate molecular weight for each compound\nd_train['MolecularWeight'] = d_train['Molecule'].apply(Descriptors.MolWt)\n\n# Print the molecular weights\nprint(d_train['MolecularWeight'])","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:54:48.969310Z","iopub.execute_input":"2023-09-16T09:54:48.970666Z","iopub.status.idle":"2023-09-16T09:54:49.219213Z","shell.execute_reply.started":"2023-09-16T09:54:48.970622Z","shell.execute_reply":"2023-09-16T09:54:49.217848Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# **Chemical Similarity And  Clustering**\n\nUse RDKit to calculate chemical similarities and a clustering algorithm to group similar compounds.","metadata":{}},{"cell_type":"code","source":"from rdkit import Chem\nfrom rdkit.Chem import DataStructs\nfrom rdkit.Chem import MACCSkeys\nfrom sklearn.cluster import AgglomerativeClustering\n\n# Convert SMILES to RDKit molecules\nd_train['Molecule'] = d_train['SMILES'].apply(Chem.MolFromSmiles)\n\n# Calculate MACCS fingerprints for each molecule\nd_train['Fingerprint'] = d_train['Molecule'].apply(MACCSkeys.GenMACCSKeys)\n\n# Calculate Tanimoto similarity between fingerprints\nsimilarity_matrix = []\nfor i, fp1 in enumerate(d_train['Fingerprint']):\n    similarities = [DataStructs.FingerprintSimilarity(fp1, fp2) for fp2 in d_train['Fingerprint']]\n    similarity_matrix.append(similarities)\n\n# Perform hierarchical clustering\nclustering = AgglomerativeClustering(n_clusters=3).fit(similarity_matrix)\nd_train['Cluster'] = clustering.labels_\n\n# Print the clustering labels\nprint(\"Clustering Labels:\")\nprint(d_train['Cluster'])","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:54:49.220789Z","iopub.execute_input":"2023-09-16T09:54:49.221241Z","iopub.status.idle":"2023-09-16T09:54:51.835864Z","shell.execute_reply.started":"2023-09-16T09:54:49.221201Z","shell.execute_reply":"2023-09-16T09:54:51.834561Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# **Enrichment Analysis**\n\nEnrichment analysis typically involves comparing a set of compounds against a database of chemical classes or substructures to identify over-represented features.","metadata":{}},{"cell_type":"code","source":"from rdkit.Chem import MACCSkeys\n\n# Define a SMARTS pattern for aromatic rings\nsmarts_pattern = Chem.MolFromSmarts('[c,C]')\n\n# Calculate MACCS fingerprints for each molecule\nd_train['Fingerprint'] = d_train['Molecule'].apply(MACCSkeys.GenMACCSKeys)\n\n# Count the number of compounds with aromatic rings\nd_train['Has_Aromatic_Ring'] = d_train['Molecule'].apply(lambda mol: mol.HasSubstructMatch(smarts_pattern))\n\n# Perform enrichment analysis\nenrichment_count = d_train['Has_Aromatic_Ring'].sum()\ntotal_compounds = len(d_train)\nenrichment_ratio = enrichment_count / total_compounds\n\nprint(f\"Enrichment of Aromatic Rings: {enrichment_ratio:.2%}\")","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:54:51.837719Z","iopub.execute_input":"2023-09-16T09:54:51.838230Z","iopub.status.idle":"2023-09-16T09:54:52.720204Z","shell.execute_reply.started":"2023-09-16T09:54:51.838180Z","shell.execute_reply":"2023-09-16T09:54:52.718276Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# **Statistical Testing**\n\nStatistical testing involves comparing groups of compounds to identify significant differences. This requires specific hypotheses and appropriate statistical tests.","metadata":{}},{"cell_type":"code","source":"from scipy.stats import ttest_ind\n\ngroup_1 = d_train[d_train['cell_type'] == 'NK cells']\ngroup_2 = d_train[d_train['cell_type'] == 'T cells CD4+']\n\nstatistic, p_value = ttest_ind(group_1['A1BG'], group_2['A1BG'])\n\nif p_value < 0.05:\n    print(f\"Statistically significant difference (p-value = {p_value:.4f})\")\nelse:\n    print(f\"No statistically significant difference (p-value = {p_value:.4f})\")\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:54:52.721952Z","iopub.execute_input":"2023-09-16T09:54:52.722479Z","iopub.status.idle":"2023-09-16T09:54:52.765148Z","shell.execute_reply.started":"2023-09-16T09:54:52.722412Z","shell.execute_reply":"2023-09-16T09:54:52.763704Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# List of columns to drop\ncolumns_to_drop = ['Molecule', 'MolecularWeight', 'Fingerprint', 'Cluster', 'Has_Aromatic_Ring']\n\n# Drop the specified columns\nd_train = d_train.drop(columns=columns_to_drop)\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:54:52.766750Z","iopub.execute_input":"2023-09-16T09:54:52.767162Z","iopub.status.idle":"2023-09-16T09:54:52.824067Z","shell.execute_reply.started":"2023-09-16T09:54:52.767126Z","shell.execute_reply":"2023-09-16T09:54:52.822794Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"d_train = pd.read_parquet('/kaggle/input/open-problems-single-cell-perturbations/de_train.parquet')","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:56:34.050018Z","iopub.execute_input":"2023-09-16T09:56:34.050541Z","iopub.status.idle":"2023-09-16T09:56:35.590144Z","shell.execute_reply.started":"2023-09-16T09:56:34.050503Z","shell.execute_reply":"2023-09-16T09:56:35.588931Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Drop duplicate rows based on all columns\nd_train.drop_duplicates(inplace=True)","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:56:39.524572Z","iopub.execute_input":"2023-09-16T09:56:39.525031Z","iopub.status.idle":"2023-09-16T09:56:43.960823Z","shell.execute_reply.started":"2023-09-16T09:56:39.524995Z","shell.execute_reply":"2023-09-16T09:56:43.959478Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# List of SMILES strings to drop\nsmiles_to_drop = [\n    'SMILES_O=C(c1ccc(Nc2nccc(-c3cc4ccccc4s3)n2)cc1)N1CCC(N2CCCC2)CC1',\n    'SMILES_O=C(c1ccc(OCCN2CCCCC2)cc1)c1c(-c2ccc(O)cc2)sc2cc(O)ccc12',\n    'SMILES_O=C1CC2(CCCC2)CC(=O)N1CCCCN1CCN(c2ncccn2)CC1',\n    'SMILES_O=C1NC(=O)C(c2cnc3ccccn23)=C1c1cn2c3c(cc(F)cc13)CN(C(=O)N1CCCCC1)CC2',\n    'SMILES_O=C1NC(=O)[C@@H](c2cn3c4c(cccc24)CCC3)[C@@H]1c1c[nH]c2ccccc12',\n    'SMILES_O=C1c2cccc3cccc(c23)C(=O)N1CCCCCC(O)=NO',\n    'SMILES_OC1(c2ccc(Cl)c(C(F)(F)F)c2)CCN(CCCC(c2ccc(F)cc2)c2ccc(F)cc2)CC1',\n    'SMILES_OCCCNc1cc(-c2ccnc(Nc3cccc(Cl)c3)n2)ccn1',\n    'SMILES_c1cc(OCCN2CCCCC2)cc(-c2[nH]nc3ccc(-c4nc[nH]n4)cc23)c1',\n    'SMILES_c1ccc2c(-c3cnn4cc(-c5ccc(N6CCNCC6)cc5)cnc34)ccnc2c1'\n]\n\n# Create a mask for rows to drop\nmask = ~d_train['SMILES'].isin(smiles_to_drop)\n\n# Apply the mask to the DataFrame\nd_train = d_train[mask]","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:57:07.130297Z","iopub.execute_input":"2023-09-16T09:57:07.130727Z","iopub.status.idle":"2023-09-16T09:57:07.178570Z","shell.execute_reply.started":"2023-09-16T09:57:07.130691Z","shell.execute_reply":"2023-09-16T09:57:07.177153Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# **Build Model and Prediction**","metadata":{}},{"cell_type":"code","source":"target_column = 'A1BG'\nX_train = d_train.drop(columns=[target_column])  # Assuming 'd_train' is your DataFrame\ny_train = d_train[target_column]","metadata":{"execution":{"iopub.status.busy":"2023-09-16T09:57:17.897999Z","iopub.execute_input":"2023-09-16T09:57:17.898487Z","iopub.status.idle":"2023-09-16T09:57:17.944499Z","shell.execute_reply.started":"2023-09-16T09:57:17.898431Z","shell.execute_reply":"2023-09-16T09:57:17.943079Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"d_train = d_train.drop(columns=['cell_type', 'sm_name'])","metadata":{"execution":{"iopub.status.busy":"2023-09-16T10:01:48.493400Z","iopub.execute_input":"2023-09-16T10:01:48.494225Z","iopub.status.idle":"2023-09-16T10:01:48.548041Z","shell.execute_reply.started":"2023-09-16T10:01:48.494180Z","shell.execute_reply":"2023-09-16T10:01:48.546672Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_train = d_train.drop('A1BG', axis=1)\ny_train = d_train['A1BG']","metadata":{"execution":{"iopub.status.busy":"2023-09-16T10:02:06.441676Z","iopub.execute_input":"2023-09-16T10:02:06.442132Z","iopub.status.idle":"2023-09-16T10:02:06.491031Z","shell.execute_reply.started":"2023-09-16T10:02:06.442094Z","shell.execute_reply":"2023-09-16T10:02:06.489577Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"categorical_columns = d_train.select_dtypes(include=['object']).columns.tolist()\n\n# Perform one-hot encoding\nd_train_encoded = pd.get_dummies(d_train, columns=categorical_columns)\n\n# Define target variable and features\nX_train = d_train_encoded.drop('A1BG', axis=1)\ny_train = d_train_encoded['A1BG']\n\n# Now, you can proceed to fit the model\nmodel = LinearRegression()\nmodel.fit(X_train, y_train)\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T10:02:13.477787Z","iopub.execute_input":"2023-09-16T10:02:13.478237Z","iopub.status.idle":"2023-09-16T10:02:17.410765Z","shell.execute_reply.started":"2023-09-16T10:02:13.478198Z","shell.execute_reply":"2023-09-16T10:02:17.409272Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.feature_selection import SelectKBest, mutual_info_regression\n\n# X_train and y_train are your training data\nnum_features_to_select = 100  # Adjust this as needed\nselector = SelectKBest(score_func=mutual_info_regression, k=num_features_to_select)\nX_train_selected = selector.fit_transform(X_train, y_train)\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T10:02:25.375807Z","iopub.execute_input":"2023-09-16T10:02:25.376405Z","iopub.status.idle":"2023-09-16T10:03:56.261779Z","shell.execute_reply.started":"2023-09-16T10:02:25.376346Z","shell.execute_reply":"2023-09-16T10:03:56.260286Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.linear_model import Ridge\n\nalpha = 1.0  # Regularization strength (adjust as needed)\nridge_model = Ridge(alpha=alpha)\nridge_model.fit(X_train_selected, y_train)\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T10:03:56.332480Z","iopub.execute_input":"2023-09-16T10:03:56.336382Z","iopub.status.idle":"2023-09-16T10:03:56.364001Z","shell.execute_reply.started":"2023-09-16T10:03:56.336245Z","shell.execute_reply":"2023-09-16T10:03:56.362372Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.decomposition import PCA\n\n#  X_train_selected from feature selection\nnum_components = 100  # Adjust this as needed\npca = PCA(n_components=num_components)\nX_train_pca = pca.fit_transform(X_train_selected)\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T10:03:56.417328Z","iopub.execute_input":"2023-09-16T10:03:56.419618Z","iopub.status.idle":"2023-09-16T10:03:56.456774Z","shell.execute_reply.started":"2023-09-16T10:03:56.419553Z","shell.execute_reply":"2023-09-16T10:03:56.455162Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.ensemble import RandomForestRegressor\n\n#  X_train_pca from dimensionality reduction\nrf_model = RandomForestRegressor(n_estimators=100, random_state=42)  # Adjust parameters as needed\nrf_model.fit(X_train_pca, y_train)\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T10:04:06.567516Z","iopub.execute_input":"2023-09-16T10:04:06.568054Z","iopub.status.idle":"2023-09-16T10:04:11.522388Z","shell.execute_reply.started":"2023-09-16T10:04:06.568006Z","shell.execute_reply":"2023-09-16T10:04:11.520916Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.model_selection import cross_val_score\n\n# Perform cross-validation\ncv_scores = cross_val_score(rf_model, X_train_pca, y_train, cv=5, scoring='neg_mean_squared_error')\n\n# Calculate RMSE\nrmse_scores = (-cv_scores)**0.5\nmean_rmse = rmse_scores.mean()\nprint(f\"Mean RMSE: {mean_rmse}\")","metadata":{"execution":{"iopub.status.busy":"2023-09-16T10:04:11.525213Z","iopub.execute_input":"2023-09-16T10:04:11.525631Z","iopub.status.idle":"2023-09-16T10:04:30.260289Z","shell.execute_reply.started":"2023-09-16T10:04:11.525593Z","shell.execute_reply":"2023-09-16T10:04:30.258947Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def mrrmse(y_true, y_pred):\n    rmse = np.sqrt(np.mean((y_true - y_pred)**2))\n    return rmse / np.mean(y_true)\n\n# The X_train_pca, y_train, and rf_model from previous steps\ncv_scores = cross_val_score(rf_model, X_train_pca, y_train, cv=5, scoring=mrrmse)\n\n# Calculate mean MRRMSE\nmean_mrrmse = np.mean(cv_scores)\nprint(f\"Mean MRRMSE: {mean_mrrmse}\")","metadata":{"execution":{"iopub.status.busy":"2023-09-16T10:04:30.261802Z","iopub.execute_input":"2023-09-16T10:04:30.262217Z","iopub.status.idle":"2023-09-16T10:04:48.874258Z","shell.execute_reply.started":"2023-09-16T10:04:30.262181Z","shell.execute_reply":"2023-09-16T10:04:48.872682Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# gene_columns is a list of gene names\ngene_columns = X_train.columns.tolist()\n\n# Create a DataFrame with zeros\nsubmission_df = pd.DataFrame(0, index=range(5), columns=['id'] + gene_columns)\n\n# Assign ids to the 'id' column\nsubmission_df['id'] = range(5)\n\n# Save the submission file\nsubmission_df.to_csv('submission.csv', index=False)\n","metadata":{"execution":{"iopub.status.busy":"2023-09-16T10:05:03.617610Z","iopub.execute_input":"2023-09-16T10:05:03.618056Z","iopub.status.idle":"2023-09-16T10:05:03.705963Z","shell.execute_reply.started":"2023-09-16T10:05:03.618019Z","shell.execute_reply":"2023-09-16T10:05:03.704543Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}