{"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":"gpu","dataSources":[{"sourceType":"competition","sourceId":67356,"databundleVersionId":8006601},{"sourceType":"datasetVersion","sourceId":8917074,"datasetId":5362625,"databundleVersionId":9078722}],"dockerImageVersionId":30615,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"#%pip install pandas==2.1.0","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:32:00.466754Z","iopub.execute_input":"2024-07-10T07:32:00.467163Z","iopub.status.idle":"2024-07-10T07:32:00.473101Z","shell.execute_reply.started":"2024-07-10T07:32:00.467129Z","shell.execute_reply":"2024-07-10T07:32:00.471748Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"This notebook can be used to train XGboost or LightGBM on different subsets of train data (with 10M subset you can run in kaggle).","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nfrom scipy import sparse\n\nfrom sklearn.model_selection import PredefinedSplit, ParameterGrid, GridSearchCV, train_test_split\nfrom sklearn.metrics import average_precision_score \n\nimport lightgbm as lgb\nimport xgboost as xgb\n\nimport json\nimport os","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:32:00.476257Z","iopub.execute_input":"2024-07-10T07:32:00.476592Z","iopub.status.idle":"2024-07-10T07:32:03.989509Z","shell.execute_reply.started":"2024-07-10T07:32:00.476565Z","shell.execute_reply":"2024-07-10T07:32:03.988368Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(pd.__version__)","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:32:03.991666Z","iopub.execute_input":"2024-07-10T07:32:03.992012Z","iopub.status.idle":"2024-07-10T07:32:03.997637Z","shell.execute_reply.started":"2024-07-10T07:32:03.991981Z","shell.execute_reply":"2024-07-10T07:32:03.996432Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Define subset name and size, features, second features, model, etc.","metadata":{}},{"cell_type":"code","source":"SUBSET_NAME = 'train_no_test_wide'\nTRAIN_SIZE = '10M'\nFEATURES = 'ecfp'\nLENGTH = '1024'\nSECOND_FEATURES = '' #mpnn_che2_10Mtnt_last15\nMODEL = 'XGb' #'lightgbm' \nPROTEINS = ['BRD4', 'HSA', 'sEH']","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:32:03.999093Z","iopub.execute_input":"2024-07-10T07:32:03.999385Z","iopub.status.idle":"2024-07-10T07:32:04.012738Z","shell.execute_reply.started":"2024-07-10T07:32:03.999359Z","shell.execute_reply":"2024-07-10T07:32:04.011651Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Create paths to input and output data (just wanted to specify variables above and have all the rest done automatically.)","metadata":{}},{"cell_type":"code","source":"DATA_PATH = '/kaggle/input/belka-data-for-gbdts'\n\nTARGETS_PATH = os.path.join(DATA_PATH, f'{SUBSET_NAME}_{TRAIN_SIZE}.parquet')\n\nTRAIN_FEATURES_FILE = f'features_{FEATURES}_{LENGTH}_{SECOND_FEATURES}_{SUBSET_NAME}_{TRAIN_SIZE}'\n\nENS_TEST_PATH = os.path.join(DATA_PATH, 'test_ensebmle_wide.parquet')\nENS_FEATURES_PATH = os.path.join(DATA_PATH, f'features_{FEATURES}_{LENGTH}_{SECOND_FEATURES}_test_ensebmle_wide')\n\nTEST_FILE = '/kaggle/input/leash-BELKA/test.csv'\nTEST_FEATURES_PATH = os.path.join(DATA_PATH,  f'features_{FEATURES}_{LENGTH}_{SECOND_FEATURES}_test')\n\nMODELS_DIR = f'{MODEL}_{FEATURES}_{LENGTH}_{SECOND_FEATURES}'\nMODELS_DIR = MODELS_DIR.replace('__', '_').rstrip('_')\nos.makedirs(MODELS_DIR, exist_ok=True)\n\nPREDS_PATH = f'{MODEL}_{FEATURES}_{LENGTH}_{SECOND_FEATURES}_{SUBSET_NAME}_{TRAIN_SIZE}_test.parquet'\nSUBMIT_PATH = f'{MODEL}_{FEATURES}_{LENGTH}_{SECOND_FEATURES}_{SUBSET_NAME}_{TRAIN_SIZE}_submit.csv'\n\ndef clean_path(*paths):\n    return [path.replace('__', '_').rstrip('_') for path in paths]\n\n(TARGETS_PATH, TRAIN_FEATURES_FILE, ENS_FEATURES_PATH, MODELS_DIR, \n PREDS_PATH, TEST_FEATURES_PATH, SUBMIT_PATH) = clean_path(\n    TARGETS_PATH, TRAIN_FEATURES_FILE, ENS_FEATURES_PATH, MODELS_DIR, \n    PREDS_PATH, TEST_FEATURES_PATH, SUBMIT_PATH\n)","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:32:04.015520Z","iopub.execute_input":"2024-07-10T07:32:04.016386Z","iopub.status.idle":"2024-07-10T07:32:04.028417Z","shell.execute_reply.started":"2024-07-10T07:32:04.016354Z","shell.execute_reply":"2024-07-10T07:32:04.027298Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Prepare the features once by loading from .npz and saving LightGBM Dataset or xgb.DMatrix into a binary file, restart the kernel.","metadata":{}},{"cell_type":"code","source":"TRAIN_FEATURES_PATH = os.path.join(DATA_PATH, f'{TRAIN_FEATURES_FILE}')\n\nif MODEL == 'lightgbm':\n    if not os.path.exists(f'{TRAIN_FEATURES_FILE}.bin'):\n        dtrain = lgb.Dataset(sparse.load_npz(f'{TRAIN_FEATURES_PATH}.npz'))  \n        dtrain.save_binary(f'{TRAIN_FEATURES_FILE}.bin')\nelif MODEL == 'XGb':\n    if not os.path.exists(f'{TRAIN_FEATURES_FILE}.buffer'):\n        dtrain = xgb.DMatrix(sparse.load_npz(f'{TRAIN_FEATURES_PATH}.npz'))\n        dtrain.save_binary(f'{TRAIN_FEATURES_FILE}.buffer')","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:32:04.029722Z","iopub.execute_input":"2024-07-10T07:32:04.030083Z","iopub.status.idle":"2024-07-10T07:32:04.043049Z","shell.execute_reply.started":"2024-07-10T07:32:04.030053Z","shell.execute_reply":"2024-07-10T07:32:04.042150Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Define a function to calculate AP scores for each protein and MAP.","metadata":{}},{"cell_type":"code","source":"def compute_ap(df, pred_cols):\n    df = pd.concat(\n        [df.loc[:, ['molecule_smiles', 'BRD4', 'HSA', 'sEH']].melt(\n            id_vars = 'molecule_smiles', \n            var_name = 'protein',\n            value_name = 'y_true'\n        ),\n         df.loc[:, ['molecule_smiles'] + pred_cols].melt(\n             id_vars = 'molecule_smiles',\n             var_name = '_',\n             value_name = 'y_pred'\n         ).drop('molecule_smiles', axis = 1)],\n        axis = 1\n    )\n    ap = df.groupby(['protein']).apply(\n        lambda x: average_precision_score(x.y_true, x.y_pred)\n    ).reset_index(name = 'AP')\n    \n    return(ap)","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:32:04.044319Z","iopub.execute_input":"2024-07-10T07:32:04.044730Z","iopub.status.idle":"2024-07-10T07:32:04.056952Z","shell.execute_reply.started":"2024-07-10T07:32:04.044700Z","shell.execute_reply":"2024-07-10T07:32:04.055902Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def compute_ap(df, pred_cols):\n    true_df = df[['molecule_smiles', 'BRD4', 'HSA', 'sEH']].melt(\n        id_vars='molecule_smiles',\n        var_name='protein',\n        value_name='y_true'\n    )\n    pred_df = df[['molecule_smiles'] + pred_cols].melt(\n        id_vars='molecule_smiles',\n        var_name='_',\n        value_name='y_pred'\n    ).drop('_', axis=1)\n    \n    combined_df = pd.concat([true_df, pred_df.drop('molecule_smiles', axis=1)], axis=1)\n    \n    ap = combined_df.groupby('protein').apply(\n        lambda x: average_precision_score(x['y_true'], x['y_pred'])\n    ).reset_index(name='AP')\n    \n    return ap","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:32:04.058326Z","iopub.execute_input":"2024-07-10T07:32:04.058608Z","iopub.status.idle":"2024-07-10T07:32:04.072913Z","shell.execute_reply.started":"2024-07-10T07:32:04.058584Z","shell.execute_reply":"2024-07-10T07:32:04.071830Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Load targets","metadata":{}},{"cell_type":"code","source":"targets = pd.read_parquet(TARGETS_PATH)[PROTEINS]\ntargets","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:32:04.074194Z","iopub.execute_input":"2024-07-10T07:32:04.074531Z","iopub.status.idle":"2024-07-10T07:32:13.908877Z","shell.execute_reply.started":"2024-07-10T07:32:04.074503Z","shell.execute_reply":"2024-07-10T07:32:13.907817Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Fit models","metadata":{}},{"cell_type":"code","source":"lgb_params = {\n    'max_depth': 11,\n    'bagging_fraction': 0.9,\n    'learning_rate': 0.05,\n    'colsample_bytree': 1,\n    'colsample_bynode': 0.5,\n    'lambda_l1': 1,\n    'objective': 'binary',\n    'lambda_l2': 1.5,\n    'num_leaves': 490,\n    'min_data_in_leaf': 50,\n    'verbose': -1,\n    'metric': 'average_precision',\n    'device': 'cpu'\n}\n\nxgb_params = {\n    'objective': 'binary:logistic',\n    'eta': 0.05,\n    'max_depth': 25,\n    'subsample': 0.2,\n    'sampling_method': 'gradient_based',\n    'colsample_bytree': 0.4,\n    'min_child_weight': 4,\n    'gamma': 2,\n    'device': 'gpu' \n}","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:32:13.910421Z","iopub.execute_input":"2024-07-10T07:32:13.910849Z","iopub.status.idle":"2024-07-10T07:32:13.918275Z","shell.execute_reply.started":"2024-07-10T07:32:13.910816Z","shell.execute_reply":"2024-07-10T07:32:13.917234Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"if MODEL == 'lightgbm':\n    dtrain = lgb.Dataset(f'{TRAIN_FEATURES_FILE}.bin') \nelif MODEL == 'XGb':\n    dtrain = xgb.DMatrix(f'{TRAIN_FEATURES_FILE}.buffer') \ndtrain","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:32:13.922307Z","iopub.execute_input":"2024-07-10T07:32:13.922771Z","iopub.status.idle":"2024-07-10T07:32:26.913847Z","shell.execute_reply.started":"2024-07-10T07:32:13.922723Z","shell.execute_reply":"2024-07-10T07:32:26.912734Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"all_indices = np.arange(len(targets))\nnp.random.seed(1)\nnp.random.shuffle(all_indices)\nvalid_idx = np.random.choice(all_indices, size = 200_000, replace = False)\ntrain_idx = np.setdiff1d(all_indices, valid_idx)\n\nprint(\"Number of samples for training\", len(train_idx))\nprint(\"Number of samples for validation:\", len(valid_idx))\nprint(\"Sanity check: intersection between train_idx and val_idx:\", np.intersect1d(train_idx, valid_idx))","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:32:26.915191Z","iopub.execute_input":"2024-07-10T07:32:26.915516Z","iopub.status.idle":"2024-07-10T07:32:29.721907Z","shell.execute_reply.started":"2024-07-10T07:32:26.915488Z","shell.execute_reply":"2024-07-10T07:32:29.720741Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"models = {}\n\nfor protein in PROTEINS:\n    \n    dtrain.set_label(targets[protein])\n\n    if MODEL == 'lightgbm':\n        bst = lgb.train(lgb_params,\n                        num_boost_round = 5000,\n                        train_set = dtrain.subset(train_idx),\n                        valid_sets = dtrain.subset(valid_idx),\n                        callbacks = [\n                            lgb.early_stopping(stopping_rounds = 30),\n                            lgb.log_evaluation(50)\n                        ]\n                       )\n        bst.save_model(os.path.join(MODELS_DIR, f'{MODEL}_model_{protein}.txt'))\n        del bst\n    elif MODEL == 'XGb':\n        evallist = [(dtrain.slice(train_idx), 'train'),\n                    (dtrain.slice(valid_idx), 'eval')]\n\n        bst = xgb.train(xgb_params,\n                        num_boost_round = 5000,\n                        dtrain = evallist[0][0],\n                        evals = evallist,\n                        early_stopping_rounds = 30,\n                        verbose_eval = 50\n                        )\n        bst.save_model(os.path.join(MODELS_DIR, f'{MODEL}_model_{protein}.ubj'))\n        del bst, evallist","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:35:19.252560Z","iopub.execute_input":"2024-07-10T07:35:19.252968Z","iopub.status.idle":"2024-07-10T07:36:15.050061Z","shell.execute_reply.started":"2024-07-10T07:35:19.252934Z","shell.execute_reply":"2024-07-10T07:36:15.048994Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Evaluate with test subset and save predictions","metadata":{}},{"cell_type":"code","source":"test_df = pd.read_parquet(ENS_TEST_PATH)\nfeatures = sparse.load_npz(f'{ENS_FEATURES_PATH}.npz')\nfeatures","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:36:35.925357Z","iopub.execute_input":"2024-07-10T07:36:35.925754Z","iopub.status.idle":"2024-07-10T07:36:43.406820Z","shell.execute_reply.started":"2024-07-10T07:36:35.925724Z","shell.execute_reply":"2024-07-10T07:36:43.405702Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"y_pred = []\nfor protein in PROTEINS:\n    if MODEL == 'lightgbm':\n        model = lgb.Booster(model_file = os.path.join(MODELS_DIR, f'{MODEL}_model_{protein}.txt'))\n        preds = model.predict(features, num_iteration = model.best_iteration)\n        y_pred.append(preds)\n        \n    elif MODEL == 'XGb':\n        model = xgb.Booster()\n        model.load_model(os.path.join(MODELS_DIR, f'{MODEL}_model_{protein}.ubj'))\n        preds = model.predict(xgb.DMatrix(features), iteration_range = (0, model.best_iteration))\n        y_pred.append(preds)","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:36:43.409045Z","iopub.execute_input":"2024-07-10T07:36:43.409879Z","iopub.status.idle":"2024-07-10T07:36:52.113356Z","shell.execute_reply.started":"2024-07-10T07:36:43.409840Z","shell.execute_reply":"2024-07-10T07:36:52.112238Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"columns_to_drop = ['buildingblock1_smiles', 'buildingblock2_smiles', 'buildingblock3_smiles', 'split2']\ntest_df = test_df.drop(columns = [col for col in columns_to_drop if col in test_df.columns])\n\npred_cols = [f'BRD4_pred_{MODEL}_{FEATURES}{LENGTH}', f'HSA_pred_{MODEL}_{FEATURES}{LENGTH}', f'sEH_pred_{MODEL}_{FEATURES}{LENGTH}']\ntest_df[pred_cols] = np.array(y_pred).T\ntest_df","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:37:04.442975Z","iopub.execute_input":"2024-07-10T07:37:04.443354Z","iopub.status.idle":"2024-07-10T07:37:04.582601Z","shell.execute_reply.started":"2024-07-10T07:37:04.443325Z","shell.execute_reply":"2024-07-10T07:37:04.581427Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import warnings\ndef filter_warnings(message, category):\n    return category is FutureWarning and \"is_sparse is deprecated and will be removed in a future version.\" in str(message)\n\nwith warnings.catch_warnings():\n    warnings.filterwarnings(\n        \"ignore\", category=FutureWarning,\n        message=\"is_sparse is deprecated and will be removed in a future version.\"\n    )\n\n    metrics = compute_ap(test_df, pred_cols)\n    metrics['index'] = 0\n\n    pivot_df = metrics.pivot(index='index', columns=['protein'], values='AP')\n    pivot_df['test_share_mean'] = metrics['AP'].mean()\n    pivot_df = pivot_df.round(4)\npivot_df","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:41:06.844250Z","iopub.execute_input":"2024-07-10T07:41:06.845114Z","iopub.status.idle":"2024-07-10T07:41:11.636688Z","shell.execute_reply.started":"2024-07-10T07:41:06.845066Z","shell.execute_reply":"2024-07-10T07:41:11.635612Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test_df.to_parquet(PREDS_PATH, index = False)","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:41:15.991888Z","iopub.execute_input":"2024-07-10T07:41:15.992283Z","iopub.status.idle":"2024-07-10T07:41:17.525451Z","shell.execute_reply.started":"2024-07-10T07:41:15.992253Z","shell.execute_reply":"2024-07-10T07:41:17.524417Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Prepare submit","metadata":{}},{"cell_type":"code","source":"test_df = pd.read_csv(TEST_FILE)\ncolumns_to_drop = ['buildingblock1_smiles', 'buildingblock2_smiles', 'buildingblock3_smiles']\ntest_df = test_df.drop(columns=[col for col in columns_to_drop if col in test_df.columns])\ntest_df","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:41:17.845160Z","iopub.execute_input":"2024-07-10T07:41:17.845574Z","iopub.status.idle":"2024-07-10T07:41:24.127266Z","shell.execute_reply.started":"2024-07-10T07:41:17.845541Z","shell.execute_reply":"2024-07-10T07:41:24.126125Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"features = sparse.load_npz(f'{TEST_FEATURES_PATH}.npz')\nfeatures","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:41:24.129553Z","iopub.execute_input":"2024-07-10T07:41:24.129952Z","iopub.status.idle":"2024-07-10T07:41:26.748778Z","shell.execute_reply.started":"2024-07-10T07:41:24.129918Z","shell.execute_reply":"2024-07-10T07:41:26.747678Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"submit = []\n\nfor protein in PROTEINS:\n    if MODEL == 'lightgbm':\n        model = lgb.Booster(model_file = os.path.join(MODELS_DIR, f'{MODEL}_model_{protein}.txt'))\n        preds = model.predict(features, num_iteration = model.best_iteration)\n        submit.append(preds)\n    elif MODEL == 'XGb':\n        model = xgb.Booster()\n        model.load_model(os.path.join(MODELS_DIR, f'{MODEL}_model_{protein}.ubj'))\n        preds = model.predict(xgb.DMatrix(features), iteration_range = (0, model.best_iteration))\n        submit.append(preds)\n        \nsubmit = np.array(submit).T\nsubmit = pd.DataFrame(submit, columns = ['BRD4', 'HSA', 'sEH'])\nsubmit = submit.reset_index(drop = True)\n\nsubmit = pd.concat([test_df['molecule_smiles'], submit], axis = 1)\nsubmit","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:41:26.750039Z","iopub.execute_input":"2024-07-10T07:41:26.750368Z","iopub.status.idle":"2024-07-10T07:41:33.937598Z","shell.execute_reply.started":"2024-07-10T07:41:26.750338Z","shell.execute_reply":"2024-07-10T07:41:33.936349Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"submit = pd.melt(\n    submit, \n    id_vars = ['molecule_smiles'], \n    value_vars = ['BRD4', 'HSA', 'sEH'], \n    value_name = 'binds', \n    var_name = 'protein_name'\n)\nsubmit = pd.merge(\n    test_df, \n    submit, \n    how = 'inner',\n    on = ['molecule_smiles', 'protein_name']\n)\nsubmit = submit[['id', 'binds']]\nsubmit = submit.drop_duplicates()\nsubmit","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:41:33.940195Z","iopub.execute_input":"2024-07-10T07:41:33.940558Z","iopub.status.idle":"2024-07-10T07:41:38.907515Z","shell.execute_reply.started":"2024-07-10T07:41:33.940527Z","shell.execute_reply":"2024-07-10T07:41:38.906275Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"submit['binds'].describe()","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:41:38.909174Z","iopub.execute_input":"2024-07-10T07:41:38.909692Z","iopub.status.idle":"2024-07-10T07:41:38.984601Z","shell.execute_reply.started":"2024-07-10T07:41:38.909624Z","shell.execute_reply":"2024-07-10T07:41:38.983428Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"submit.to_csv(SUBMIT_PATH, index=False)","metadata":{"execution":{"iopub.status.busy":"2024-07-10T07:41:38.986312Z","iopub.execute_input":"2024-07-10T07:41:38.986757Z","iopub.status.idle":"2024-07-10T07:41:44.935048Z","shell.execute_reply.started":"2024-07-10T07:41:38.986720Z","shell.execute_reply":"2024-07-10T07:41:44.933870Z"},"trusted":true},"execution_count":null,"outputs":[]}]}