{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.13","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":9012943,"sourceType":"datasetVersion","datasetId":5430523},{"sourceId":9031752,"sourceType":"datasetVersion","datasetId":5443798},{"sourceId":9032701,"sourceType":"datasetVersion","datasetId":5444531}],"dockerImageVersionId":30698,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Installing and importing ","metadata":{}},{"cell_type":"code","source":"#install libraries\n\n!pip install mapply\n!pip install rdkit\n!pip install duckdb","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#import libraries\n\nimport dask\nimport dask.dataframe as dd\nimport os\nimport pickle\nimport duckdb\nimport mapply\nfrom joblib import Parallel, delayed\nimport time\nfrom tqdm import tqdm\nimport pyarrow.parquet as pq\n\n#enable parallel execution\nmapply.init(n_workers=-1, progressbar=True)\n\nimport pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nimport seaborn as sns\n\nfrom rdkit import Chem\nfrom rdkit.Chem import Draw\nfrom rdkit.Chem import AllChem\n# Suppress RDKit warnings including deprecation warnings\nfrom rdkit import RDLogger\nRDLogger.DisableLog('rdApp.*')\n\nfrom catboost import CatBoostClassifier\nfrom sklearn.metrics import confusion_matrix, ConfusionMatrixDisplay, accuracy_score, precision_score, f1_score\nfrom sklearn.model_selection import train_test_split, cross_val_score\nfrom sklearn.metrics import make_scorer, average_precision_score\nfrom sklearn.preprocessing import OneHotEncoder\n\nimport warnings\nwarnings.filterwarnings('ignore')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Splitting Into Training And Validation Sets ","metadata":{}},{"cell_type":"code","source":"df_train_path = '/kaggle/input/leash-BELKA/train.parquet'\nsplit_cv_path = '/kaggle/input/belka-cv-split/train_folds.parquet'\n\n# Read the Parquet file\ndf_train = dd.read_parquet(df_train_path)\ndf_cv = dd.read_parquet(split_cv_path)\n\ndisplay(df_train.head())\ndisplay(df_cv.head())","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Add binds column\ndf_cv['binds'] = (df_cv['binds_BRD4'] == 1) | (df_cv['binds_HSA'] == 1) | (df_cv['binds_sEH'] == 1)\ndf_cv['binds'] = df_cv['binds'].astype(int)\n\ndf_cv.head()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Filter the DataFrame for false values in 'full_test' column (training)\nfiltered_train_df = df_cv[df_cv['full_test'] == False]\n\n# Filter the DataFrame for true values in 'full_test' column (validation)\nfiltered_val_df = df_cv[df_cv['full_test'] == True]\n\ndisplay(filtered_train_df.head())\ndisplay(filtered_val_df.head())","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Sampling And Preprocessing","metadata":{}},{"cell_type":"code","source":"filtered_train_df.compute().to_parquet('train_final.parquet', engine='pyarrow')","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_path = \"/kaggle/working/train_final.parquet\"\n\n# Create a DuckDB in-memory connection\ncon = duckdb.connect()\n\n# Perform the query to sample the data\ndf_sampled = con.query(f\"\"\"\n    (SELECT *\n     FROM parquet_scan('{train_path}')\n     WHERE binds = 0\n     ORDER BY RANDOM()\n     LIMIT 835000)\n    UNION ALL\n    (SELECT *\n     FROM parquet_scan('{train_path}')\n     WHERE binds = 1\n     ORDER BY RANDOM()\n     LIMIT 835000)\n\"\"\").df()\n\n# Close the connection\ncon.close()\n\n# Shuffle the final sampled DataFrame\ndf_final = df_sampled.sample(frac=1, random_state=42).reset_index(drop=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Perform a merge to filter based on the molecule_smiles in df_sampled\ndf_training = df_train.merge(df_final[['molecule_smiles']], on='molecule_smiles', how='inner')\ndf = df_training.compute()\ndf.to_parquet('train_5M.parquet', engine='pyarrow')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Function to process each chunk\ndef process_chunk(chunk):\n    # Convert SMILES to RDKit molecules\n    chunk['molecule'] = chunk['molecule_smiles'].apply(Chem.MolFromSmiles)\n    \n    # Generate ECFPs\n    def generate_ecfp(molecule, radius=3, bits=256):\n        if molecule is None:\n            return None\n        return list(AllChem.GetMorganFingerprintAsBitVect(molecule, radius, nBits=bits))\n    \n    chunk['ecfp'] = chunk['molecule'].apply(generate_ecfp)\n    \n    return chunk['ecfp']\n\n# Function to process the DataFrame in chunks with a progress bar\ndef process_dataframe_in_chunks(df, chunksize=100000):\n    num_chunks = (len(df) + chunksize - 1) // chunksize  # Calculate number of chunks\n\n    for i in tqdm(range(num_chunks), desc=\"Processing chunks\", unit=\"chunk\"):\n        start = i * chunksize\n        end = min((i + 1) * chunksize, len(df))\n        chunk = df.iloc[start:end].copy()  # Use .copy() to avoid SettingWithCopyWarning\n        df.iloc[start:end, df.columns.get_loc('ecfp')] = process_chunk(chunk)\n\n\n# Initialize the 'ecfp' column\ndf['ecfp'] = None\n\n# Measure the total time taken for processing\nstart_time = time.time()\n\n# Process the DataFrame in chunks\nprocess_dataframe_in_chunks(df)\n\n# Calculate the total time taken\nend_time = time.time()\ntotal_time = end_time - start_time\n\n# Print the total time taken\nprint(f\"Total time taken: {total_time:.2f} seconds\")","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Training","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nfrom sklearn.preprocessing import OneHotEncoder\nfrom catboost import CatBoostClassifier\nimport pickle\nfrom tqdm import tqdm\nfrom sklearn.metrics import average_precision_score\n\n# Assuming df is your processed DataFrame\nchunksize = 100000  # Adjust based on memory capacity\n\n# Initialize the one-hot encoder\nonehot_encoder = OneHotEncoder(sparse=False)\nonehot_encoder.fit(np.array(df['protein_name']).reshape(-1, 1))\n\n# Save the one-hot encoder\nwith open('onehot_encoder_10M.pkl', 'wb') as file:\n    pickle.dump(onehot_encoder, file)\n    \n# Initialize the CatBoostClassifier\ncatboost_model = CatBoostClassifier(random_state=42, verbose=0)\n\n# Function to process each chunk\ndef process_chunk(chunk):\n    # One-hot encode the 'protein_name' column\n    protein_onehot = onehot_encoder.transform(np.array(chunk['protein_name']).reshape(-1, 1))\n    \n    # Combine ECFPs and one-hot encoded 'protein_name'\n    X_chunk = [np.concatenate((ecfp, protein)) for ecfp, protein in zip(chunk['ecfp'].tolist(), protein_onehot.tolist())]\n    y_chunk = chunk['binds'].tolist()\n    \n    return X_chunk, y_chunk\n\n# Process and train on all data in chunks with a progress bar\nnum_chunks = (len(df) + chunksize - 1) // chunksize  # Calculate number of chunks\nfor start in tqdm(range(0, len(df), chunksize), total=num_chunks, desc=\"Training CatBoost Model\", unit=\"chunk\"):\n    end = start + chunksize\n    df_chunk = df.iloc[start:end]\n    X_chunk, y_chunk = process_chunk(df_chunk)\n    catboost_model.fit(X_chunk, y_chunk, init_model=catboost_model if start > 0 else None)\n\n# Define path to save the fully trained model\nmodel_path_all = 'catboost_model_all_10M.pkl'\n\n# Save the fully trained model using pickle\nwith open(model_path_all, 'wb') as f:\n    pickle.dump(catboost_model, f)\n\nprint(f\"Fully trained CatBoost model saved to: {model_path_all}\")    \n\n# Make predictions on the entire dataset in chunks\ny_true_all = []\ny_pred_proba_all = []\n\nfor start in tqdm(range(0, len(df), chunksize), total=num_chunks, desc=\"Making Predictions\", unit=\"chunk\"):\n    end = start + chunksize\n    df_chunk = df.iloc[start:end]\n    X_chunk, y_chunk = process_chunk(df_chunk)\n    y_pred_proba_chunk = catboost_model.predict_proba(X_chunk)[:, 1]\n    y_true_all.extend(y_chunk)\n    y_pred_proba_all.extend(y_pred_proba_chunk)\n\n# Calculate the mean average precision on the entire dataset\nall_map_score = average_precision_score(y_true_all, y_pred_proba_all)\nprint(f\"CatBoost Mean Average Precision (mAP) on Entire Data: {all_map_score:.2f}\")\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Validation Dataframe","metadata":{}},{"cell_type":"code","source":"filtered_val_df.compute().to_parquet('val_final.parquet', engine='pyarrow')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"val_path = \"/kaggle/working/val_final.parquet\"\n\n# Create a DuckDB in-memory connection\ncon = duckdb.connect()\n\n# Perform the query to sample the data\ndf_sampled = con.query(f\"\"\"\n    (SELECT *\n     FROM parquet_scan('{val_path}')\n     WHERE binds = 0\n     ORDER BY RANDOM()\n     LIMIT 2000000)\n    UNION ALL\n    (SELECT *\n     FROM parquet_scan('{val_path}')\n     WHERE binds = 1\n     ORDER BY RANDOM()\n     LIMIT 1000000)\n\"\"\").df()\n\n# Close the connection\ncon.close()\n\n# Shuffle the final sampled DataFrame\ndf_final_val = df_sampled.sample(frac=1, random_state=42).reset_index(drop=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Perform a merge to filter based on the molecule_smiles in df_sampled\ndf_validation = df_train.merge(df_final_val[['molecule_smiles']], on='molecule_smiles', how='inner')\ndf_validation.compute().to_parquet('final_val.parquet', engine='pyarrow')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Inference","metadata":{}},{"cell_type":"code","source":"chunksize = 100000 \n\n# Specify the path and filename for the test file\ntest_file = '/kaggle/input/validation/final_val (1).parquet'\n\n# Specify the path and filename for the model file\nmodel_file = '/kaggle/input/cat-1m/catboost_model_all.pkl'\n\n# Load the pre-trained CatBoost model from a pickle file\nwith open(model_file, 'rb') as file:\n    catboost_model = pickle.load(file)\n\n# Function to generate ECFPs\ndef generate_ecfp(molecule, radius=3, bits=256):\n    if molecule is None:\n        return None\n    return list(AllChem.GetMorganFingerprintAsBitVect(molecule, radius, nBits=bits))\n\n# Load the one-hot encoder\nwith open('/kaggle/input/cat-1m/onehot_encoder.pkl', 'rb') as file:\n    onehot_encoder = pickle.load(file)\n\n# Read the entire parquet file\ndf_test_full = pd.read_parquet(test_file, engine=\"pyarrow\")\n\n# Process the test data chunk by chunk with a progress bar\nnum_chunks = (len(df_test_full) + chunksize - 1) // chunksize  # Calculate number of chunks\n\nall_probabilities = []\nall_true_labels = []\n\nfor start in tqdm(range(0, len(df_test_full), chunksize), total=num_chunks, desc=\"Processing chunks\", unit=\"chunk\"):\n    end = start + chunksize\n    df_test = df_test_full.iloc[start:end]\n    \n    # Generate ECFPs for the molecule_smiles\n    df_test['molecule'] = df_test['molecule_smiles'].apply(Chem.MolFromSmiles)\n    df_test['ecfp'] = df_test['molecule'].apply(generate_ecfp)\n\n    # One-hot encode the protein_name\n    protein_onehot = onehot_encoder.transform(df_test['protein_name'].values.reshape(-1, 1))\n\n    # Combine ECFPs and one-hot encoded protein_name\n    X_test = [ecfp + protein.tolist() for ecfp, protein in zip(df_test['ecfp'].tolist(), protein_onehot)]\n\n    # Predict the probabilities\n    probabilities = catboost_model.predict_proba(X_test)[:, 1]\n    \n    # Collect probabilities and true labels\n    all_probabilities.extend(probabilities)\n    all_true_labels.extend(df_test['binds'].tolist())  # Assuming 'binds' column exists in the test data","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Results","metadata":{}},{"cell_type":"code","source":"# Calculate the mean average precision on the entire test dataset\ninference_map_score = average_precision_score(all_true_labels, all_probabilities)\nprint(f\"CatBoost Mean Average Precision (mAP) on Test Data: {inference_map_score:.2f}\")","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nimport numpy as np\n\n# Data\ncategories = ['1 Million', '5 Million', '60 Million\\n(ensembled)']\nvalues1 = [0.72, 0.63, 0.85]\nvalues2 = [0.43, 0.38, 0.52]\n\n\n# Number of categories\nn = len(categories)\n\n# Creating the bar width\nbar_width = 0.35\n\n# X positions for the bars\nr1 = np.arange(n)\nr2 = [x + bar_width for x in r1]\n\n# Plotting\nfig, ax = plt.subplots()\nbars1 = ax.bar(r1, values1, color='darkblue', width=bar_width, label='Training')\nbars2 = ax.bar(r2, values2, color='blue', width=bar_width, label='Validation')\n\n# Adding labels\nax.set_xlabel('Trained Data')\nax.set_ylabel('MAP Score')\nax.set_title('Model variations')\nax.set_xticks([r + bar_width / 2 for r in range(n)])\nax.set_xticklabels(categories)\n\n# Adding value labels on the bars\nfor bar in bars1:\n    yval = bar.get_height()\n    ax.text(bar.get_x() + bar.get_width() / 2, yval, yval, ha='center', va='bottom', fontsize=10)\nfor bar in bars2:\n    yval = bar.get_height()\n    ax.text(bar.get_x() + bar.get_width() / 2, yval, yval, ha='center', va='bottom', fontsize=10)\n\n# Adding legend\nax.legend()\n\n# Removing top and right borders\nax.spines['top'].set_visible(False)\nax.spines['right'].set_visible(False)\n\n# Show the plot\nplt.show()\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Calculate the confusion matrix\nthreshold = 0.5  # You can adjust the threshold based on your requirement\npredicted_labels = [1 if prob >= threshold else 0 for prob in all_probabilities]\nconf_matrix = confusion_matrix(all_true_labels, predicted_labels)\n\n# Normalize the confusion matrix by row (true labels)\nconf_matrix_normalized = conf_matrix.astype('float') / conf_matrix.sum(axis=1)[:, np.newaxis]\n\n# Display the normalized confusion matrix with percentages\ndisp = ConfusionMatrixDisplay(confusion_matrix=conf_matrix_normalized, display_labels=[\"Negative\", \"Positive\"])\ndisp.plot(cmap=plt.cm.Blues, values_format=\".2f\")\nplt.title(\"Confusion Matrix_1M\")\nplt.show()\n\n# Calculate additional metrics\naccuracy = accuracy_score(all_true_labels, predicted_labels)\nprecision = precision_score(all_true_labels, predicted_labels)\nf1 = f1_score(all_true_labels, predicted_labels)\n\nprint(f\"Accuracy: {accuracy:.2f}\")\nprint(f\"Precision: {precision:.2f}\")\nprint(f\"F1 Score: {f1:.2f}\")","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}