{"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":67356,"databundleVersionId":8006601,"sourceType":"competition"},{"sourceId":176652274,"sourceType":"kernelVersion"}],"dockerImageVersionId":30684,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Leash Bio - Predict New Medicines with BELKA\n\nAfter 4 weeks, I'm still stuck with the CV strategy.\nI created this Notebook to test different CV strategies.\n1. Sample data. All positve cases plus 8 millions rows for no positives(9,5M rows * 3 Targets).\n2. Simple Xgb model, no tunning and no GPU.\n3. Different model for each target..\n4. Split in two different problems, shared data (tests smiles exist in train) and nonshared data (tests smiles no exist in train).\n  1. For shared data, suffle rows and split. 9M rows for trainig 500K for validation.\n  2. For nonshare data,four folds of unique smiles, non shared bettein train and validation and rows with at least one positive case.\n  \n  \nVersion:\n* V.8 - VarianceThreshold. select columns with Variance Threshold  0.01 ->training with 940 columns.\n* V.9 - New down Sample data set, 4,5M rows(created like previous one of 9,5M). XGB model training iterations from 3000 to 1000. Variance Threshold 947 rows.\n* V.11 - Back to 9,5M rows and 3000 iterations for training.\n* ....\n* ... some experiments\n* ....\n* V.22 - VarianceThreshold of 0.005 and training with 1373 features. One XGB model for each target training iterations to 4000 and adding 'colsample_bytree':0.8'.\n\n* V.23 - VarianceThreshold of 0.005 and training with 1373 features. One XGB model for each target training iterations to 4000 and adding 'colsample_bytree':0.6'.","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2024-04-29T07:47:58.166649Z","iopub.execute_input":"2024-04-29T07:47:58.167386Z","iopub.status.idle":"2024-04-29T07:48:00.866609Z","shell.execute_reply.started":"2024-04-29T07:47:58.167352Z","shell.execute_reply":"2024-04-29T07:48:00.865726Z"}}},{"cell_type":"code","source":"import numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport pickle\nimport random, os, gc\nfrom scipy import sparse\nfrom sklearn.metrics import average_precision_score\n\nfrom xgboost import DMatrix\nimport xgboost as xgb\n\n\n# Sample data based on https://www.kaggle.com/datasets/shlomoron/belka-shrunken-train-set","metadata":{"execution":{"iopub.status.busy":"2024-05-31T06:04:59.335345Z","iopub.execute_input":"2024-05-31T06:04:59.336412Z","iopub.status.idle":"2024-05-31T06:05:02.314511Z","shell.execute_reply.started":"2024-05-31T06:04:59.336372Z","shell.execute_reply":"2024-05-31T06:05:02.313540Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def seed_everything(seed: int):\n    random.seed(seed)\n    np.random.seed(seed)\n    os.environ['PYTHONHASHSEED'] = str(seed) \n    \nseed_everything(42)","metadata":{"execution":{"iopub.status.busy":"2024-05-31T06:05:02.316461Z","iopub.execute_input":"2024-05-31T06:05:02.317346Z","iopub.status.idle":"2024-05-31T06:05:02.325524Z","shell.execute_reply.started":"2024-05-31T06:05:02.317305Z","shell.execute_reply":"2024-05-31T06:05:02.323869Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def fold_properties(df):\n    print('Shape :',df.shape)\n    print('block1 :',df.buildingblock1_smiles.nunique())\n    print('block2 :',df.buildingblock2_smiles.nunique())\n    print('block3 :',df.buildingblock3_smiles.nunique())    \n    for target in TARGETS:\n        print(f'Positive Cases Ratio for {target} :', target,(df[target].sum()/df.shape[0]))","metadata":{"execution":{"iopub.status.busy":"2024-05-31T06:05:02.331857Z","iopub.execute_input":"2024-05-31T06:05:02.332278Z","iopub.status.idle":"2024-05-31T06:05:02.338921Z","shell.execute_reply.started":"2024-05-31T06:05:02.332246Z","shell.execute_reply":"2024-05-31T06:05:02.338041Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"TARGETS=['binds_BRD4','binds_HSA','binds_sEH']\ndtypes = {'buildingblock1_smiles': np.int16, 'buildingblock2_smiles': np.int16, 'buildingblock3_smiles': np.int16,\n          'binds_BRD4':np.byte, 'binds_HSA':np.byte, 'binds_sEH':np.byte}\n\ntrain = pd.read_csv('/kaggle/input/leashbio-10m-data-sample/train_sample.csv', dtype = dtypes, usecols=[0,1,2,4,5,6] )\ntrain_index=train.index.tolist()\nrandom.shuffle(train_index)\nprint('Train Sample Properties')\nfold_properties(train)\n    \n#Load train_ecfp sparce matrix\ntrain_ecfp = sparse.load_npz(\"/kaggle/input/leashbio-10m-data-sample/train_ecfp.npz\")\nprint('train_ecfp Shape :',train_ecfp.shape)\n\ntrain_B1_unique=train.buildingblock1_smiles.unique().tolist()\ntrain_B2_unique=train.buildingblock2_smiles.unique().tolist()\ntrain_B3_unique=train.buildingblock3_smiles.unique().tolist()\ntrain_unique_smiles=set(train_B1_unique+train_B2_unique+train_B3_unique)\n","metadata":{"execution":{"iopub.status.busy":"2024-05-31T06:05:02.340436Z","iopub.execute_input":"2024-05-31T06:05:02.341293Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.feature_selection import VarianceThreshold\nfor i in range(6):\n    threshold=0.1/10**i\n    var_thresh = VarianceThreshold(threshold=threshold)\n    var_thresh.fit(train_ecfp[:100000].A)\n    var_thresh_index=var_thresh.get_support()\n    print(threshold,sum(var_thresh_index))\n\n\n#select columns with Variance Threshold  001]\n    \n#Apply mask to ecfp sparse matrix.\nvar_thresh = VarianceThreshold(threshold=0.005)\nvar_thresh.fit(train_ecfp[:100000].A)\nvar_thresh_index_1=var_thresh.get_support()\ntrain_ecfp_1=train_ecfp[:,var_thresh_index_1]\nprint('train_ecfp Shape :',train_ecfp_1.shape)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Shared smiles Models","metadata":{"execution":{"iopub.status.busy":"2024-04-29T08:06:00.927776Z","iopub.execute_input":"2024-04-29T08:06:00.928198Z","iopub.status.idle":"2024-04-29T08:06:00.932687Z","shell.execute_reply.started":"2024-04-29T08:06:00.928162Z","shell.execute_reply":"2024-04-29T08:06:00.931552Z"}}},{"cell_type":"code","source":"iterations = 4000\nearly_stopping_rounds = 100\nverbose_eval = 100\n\nxgb_models=[]\nvalid_preds_score=[]\n    \nprint('#*****#')\nprint(f'train Data')\ntrain_fold_index=train_index[:9*10**6]\nfold_properties(train.loc[train_fold_index])\nX_tr=train_ecfp_1[(train_fold_index),:]\n\n\nprint('#*****#')\nprint('Valid Data')\nvalid_fold_index=train_index[9*10**6:]\nfold_properties(train.loc[valid_fold_index])\nX_val=train_ecfp_1[(valid_fold_index),:]\n    \nfor i,target in enumerate(TARGETS):\n    print(f'training target {target}')\n    y_tr=train.loc[train_fold_index][target].values\n    y_val=train.loc[valid_fold_index][target].values \n        \n    xtrain = DMatrix(data=X_tr, label=y_tr)\n    xvalid = DMatrix(data=X_val, label=y_val)\n        \n    scale_pos_weight=(len(y_tr)-y_tr.sum())/y_tr.sum()\n    print('scale_pos_weight :',scale_pos_weight)\n    \n    xgb_params= { 'objective':'binary:logistic', 'eval_metric':'aucpr',\n                 'learning_rate': 0.2, \"n_jobs\":12,\"seed\": 42,'colsample_bytree':0.6, 'scale_pos_weight':scale_pos_weight}\n\n    XModel = xgb.train(xgb_params, xtrain,\n                           evals=[(xvalid,'validation')],\n                           verbose_eval=verbose_eval,\n                           early_stopping_rounds=early_stopping_rounds,\n                           xgb_model=None,\n                           num_boost_round=iterations)\n        \n\n    y_pred_proba = XModel.predict(xvalid)  # Probability of the positive class\n    map_score = average_precision_score(y_val, y_pred_proba)\n    valid_preds_score.append(map_score)\n    print(f\"Mean Average Precision for Valid {target}: {map_score:.2f}\")\n    \n    xgb_models.append(XModel)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('binds_BRD4 :',np.mean(valid_preds_score[0::3]))\nprint('binds_HSA :' ,np.mean(valid_preds_score[1::3]))\nprint('binds_sEH :' ,np.mean(valid_preds_score[2::3]))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"del train,xtrain,X_tr,xvalid,X_val,y_tr,y_val,valid_preds_score\ngc.collect()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"preds_cols=['BRD4','HSA','sEH']\ntest=pd.read_parquet('/kaggle/input/leash-BELKA/test.parquet', engine = 'pyarrow')\nblocks_dict=np.load('/kaggle/input/leashbio-10m-data-sample/blocks_dict.npy', allow_pickle=True)\nblocks_dict = blocks_dict.tolist()\ntest['buildingblock1_smiles']=test['buildingblock1_smiles'].map(blocks_dict).values.astype('uint16')\ntest['buildingblock2_smiles']=test['buildingblock2_smiles'].map(blocks_dict).values.astype('uint16')\ntest['buildingblock3_smiles']=test['buildingblock3_smiles'].map(blocks_dict).values.astype('uint16')\n\ntest_preds=pd.DataFrame(test['molecule_smiles'].unique(),columns=['molecule_smiles'])\ntest_preds=test[['molecule_smiles','buildingblock1_smiles','buildingblock2_smiles','buildingblock3_smiles']].copy()\ntest_preds.drop_duplicates(inplace=True)\ntest_preds=test_preds.reset_index(drop=True) \n\ntest_ecfp=sparse.load_npz(\"/kaggle/input/leashbio-10m-data-sample/test_ecfp.npz\")\ntest_ecfp_1=test_ecfp[:,var_thresh_index_1]\nprint('test_ecfp Shape :',train_ecfp_1.shape)\n\n\ntest_B1_unique=test_preds.buildingblock1_smiles.unique().tolist()\ntest_B2_unique=test_preds.buildingblock2_smiles.unique().tolist()\ntest_B3_unique=test_preds.buildingblock3_smiles.unique().tolist()\ntest_unique_smiles=set(test_B1_unique+test_B2_unique+test_B3_unique)\n\ntest_preds[preds_cols]=np.zeros((len(test_preds),3))\n\nprint('Shape Test :', test.shape)\ntest_no_overlap=[bb for bb in test_unique_smiles if bb not in train_unique_smiles]\ntrain_no_overlap=[bb for bb in train_unique_smiles if bb not in test_unique_smiles]\nprint('Test Block Unique :', len(test_unique_smiles))\nprint('Train Block Unique :', len(train_unique_smiles))\nprint('Test Block not in train :', len(test_no_overlap))\nprint('Train Block not in test :', len(train_no_overlap))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for i,target in enumerate(TARGETS):    \n    test_target=target.split('_')[1]      \n    test_preds[test_target]=xgb_models[i].predict(DMatrix(data=test_ecfp_1))\n\ntest_BRD4=test_preds[['molecule_smiles','BRD4']].copy()\ntest_BRD4['protein_name']='BRD4'\ntest_BRD4=test_BRD4.rename(columns={\"BRD4\": \"binds_1\"})\n\ntest_HSA=test_preds[['molecule_smiles','HSA']].copy()\ntest_HSA['protein_name']='HSA'\ntest_HSA=test_HSA.rename(columns={\"HSA\": \"binds_1\"})\n\ntest_sEH=test_preds[['molecule_smiles','sEH']].copy()\ntest_sEH['protein_name']='sEH'\ntest_sEH=test_sEH.rename(columns={\"sEH\": \"binds_1\"})\n\ntest_preds_1=pd.concat([test_BRD4,test_HSA,test_sEH])\ntest=pd.merge(test,test_preds_1, on=['molecule_smiles','protein_name'],how = 'left')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sample_submission = pd.read_csv(\"/kaggle/input/leash-BELKA/sample_submission.csv\")\nsample_submission['binds']=test['binds_1']\nsample_submission.to_csv('submission.csv', index=False)\nsample_submission.head()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(sample_submission.binds.max())\nprint(sample_submission.binds.min())\nprint(sample_submission.binds.mean())","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}