{"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":"gpu","dataSources":[{"sourceId":67356,"databundleVersionId":8006601,"sourceType":"competition"},{"sourceId":8519891,"sourceType":"datasetVersion","datasetId":4784530}],"dockerImageVersionId":30734,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"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":"2024-06-30T23:09:50.312553Z","iopub.execute_input":"2024-06-30T23:09:50.313340Z","iopub.status.idle":"2024-06-30T23:09:50.682443Z","shell.execute_reply.started":"2024-06-30T23:09:50.313201Z","shell.execute_reply":"2024-06-30T23:09:50.681523Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Leash - ECFPs and 1dCNN","metadata":{}},{"cell_type":"markdown","source":"Inspired by [AH](https://www.kaggle.com/code/ahmedelfazouan/belka-1dcnn-starter-with-all-data/notebook)\nand [ANDREW D. BLEVINS](https://www.kaggle.com/code/andrewdblevins/leash-tutorial-ecfps-and-random-forest)\n\nIdea:\n* Perform the simplist conversion from SMILES to ECFPs [1](https://pubs.acs.org/doi/10.1021/ci100050t) with 2048 bits. \n    * Due to hash collisions, multiple molecules can theoretically generate the same fingerprint. However, Morgan fingerprints with 2048 bits are quite sparse, reducing the likelihood of such collisions.\n* Train a 1dcnn model on 20 epochs (5 epochs for showcasing).\n* All variable are modify for easier version-saving.\n* Run-time for all data on TPU with preprocessed data is around 4-5hrs (on Colab).\n    * Have trouble running TPU on Kaggle.\n    * Warning: data Processing may take up to hours.    \n\n\n**When training with all data, the LB score comes up to 0.32.**\nThis notebook does not produce correct output!","metadata":{}},{"cell_type":"markdown","source":"## Environment Configuration","metadata":{}},{"cell_type":"code","source":"import sys\nimport random\nimport tensorflow as tf\n\nprint('Tensorflow version ' + tf.__version__)","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:09:50.683778Z","iopub.execute_input":"2024-06-30T23:09:50.684522Z","iopub.status.idle":"2024-06-30T23:09:53.888365Z","shell.execute_reply.started":"2024-06-30T23:09:50.684486Z","shell.execute_reply":"2024-06-30T23:09:53.887390Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class CFG:\n    PREPROCESS = True\n    EPOCHS = 5 #20\n    BATCH_SIZE = 256 #4096\n    LR = 1e-3\n    WD = 0.05\n\n    NBR_FOLDS = 15\n    SELECTED_FOLDS = [0]\n\n    SEED = 2024\n    INP_LEN = 256","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:09:53.889566Z","iopub.execute_input":"2024-06-30T23:09:53.890082Z","iopub.status.idle":"2024-06-30T23:09:53.895279Z","shell.execute_reply.started":"2024-06-30T23:09:53.890054Z","shell.execute_reply":"2024-06-30T23:09:53.894278Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Detect hardware, return appropriate distribution strategy\ntry:\n    # connect to a TPU and instantiate a distribution strategy\n    tpu = tf.distribute.cluster_resolver.TPUClusterResolver(tpu='local')\n    tf.config.experimental_connect_to_cluster(tpu)\n    tf.tpu.experimental.initialize_tpu_system(tpu)\n    strategy = tf.distribute.TPUStrategy(tpu)\n    all_data = True\n    \n    print(\"Running on TPU\")\n    print(\"REPLICAS: \", tpu_strategy.num_replicas_in_sync)\nexcept tf.errors.NotFoundError:\n    all_data = False\n    strategy = tf.distribute.OneDeviceStrategy(device=\"/gpu:0\")\n    print(\"Not on TPU\")","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:09:53.899385Z","iopub.execute_input":"2024-06-30T23:09:53.899669Z","iopub.status.idle":"2024-06-30T23:09:54.230374Z","shell.execute_reply.started":"2024-06-30T23:09:53.899645Z","shell.execute_reply":"2024-06-30T23:09:54.229422Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def set_seeds(seed):\n    os.environ['PYTHONHASHSEED'] = str(seed)\n    random.seed(seed)\n    tf.random.set_seed(seed)\n    np.random.seed(seed)\n\nset_seeds(seed=CFG.SEED)","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:09:54.231596Z","iopub.execute_input":"2024-06-30T23:09:54.231966Z","iopub.status.idle":"2024-06-30T23:09:54.237876Z","shell.execute_reply.started":"2024-06-30T23:09:54.231931Z","shell.execute_reply":"2024-06-30T23:09:54.236977Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Data Preprocessing","metadata":{}},{"cell_type":"markdown","source":"### Data Minimization and SMILE Conversion (ECFP Approach)","metadata":{}},{"cell_type":"code","source":"if CFG.PREPROCESS:\n    !pip install fastparquet -q\n    !pip install rdkit\n    \n    from rdkit import Chem\n    from rdkit.Chem import AllChem\n\n    # Generate ECFPs\n    def generate_ecfp(smiles, radius=2, bits=CFG.INP_LEN*8):\n        ecfp_arr = []\n        sm_count = 0\n        for s in smiles:\n            # Convert SMILES to RDKit molecules\n            mol = Chem.MolFromSmiles(s)\n            # Convert RDKit molecules to ECFP\n            fp_gen = Chem.rdFingerprintGenerator.GetMorganGenerator(radius=radius, fpSize=bits)\n    #         ecfp_lst.append(list(AllChem.GetMorganFingerprintAsBitVect(mol, radius, nBits=bits)))\n            ecfp = list(fp_gen.GetFingerprint(mol))\n            # Repack bits\n            ecfp_arr.append(np.packbits(ecfp))\n            \n            sm_count += 1\n            done = int(50 * sm_count / int(len(smiles)))\n            sys.stdout.write(f\"\\r[{'=' * done}{' ' * (50-done)}] {sm_count} smiles processed\")\n            sys.stdout.flush()\n\n        return ecfp_arr\n\n    if all_data:\n        train_path = '/kaggle/input/leash-BELKA/train.parquet'\n        train_raw = pd.read_parquet(train_path)\n    else:\n        train_path = '/kaggle/input/leash-BELKA/train.csv'\n        train_raw = pd.read_csv(train_path, nrows=90000)\n    \n    protein = train_raw['protein_name'].unique()\n\n    # Every smile comes in a set of three (protein binds)\n    # Only process the unique smiles\n    train_smiles = train_raw[train_raw['protein_name']==protein[0]]['molecule_smiles'].values\n    assert (train_smiles!=train_raw[train_raw['protein_name']==protein[1]]['molecule_smiles'].values).sum() == 0\n    assert (train_smiles!=train_raw[train_raw['protein_name']==protein[2]]['molecule_smiles'].values).sum() == 0\n\n#     # For preprocessing (comment out due to excessive run-time)\n#     test_path = '/kaggle/input/leash-BELKA/test.parquet'\n#     test_raw = pd.read_parquet(test_path)\n#     test_smiles = test_raw['molecule_smiles'].unique()\n\n    # ECFP lists and protein binds files\n    train_ecfp = generate_ecfp(train_smiles)\n#     test_ecfp = generate_ecfp(test_smiles)\n#     np.savez('train.ecfp4.packed.npz', ecfp=train_ecfp)\n#     np.savez('test.ecfp4.packed.npz', ecfp=test_ecfp)\n\n    bind1 = train_raw[train_raw['protein_name']==protein[0]]['binds'].values\n    bind2 = train_raw[train_raw['protein_name']==protein[1]]['binds'].values\n    bind3 = train_raw[train_raw['protein_name']==protein[2]]['binds'].values\n    train_bind = list(zip(bind1, bind2, bind3))\n#     np.savez('train.bind.npz', bind=train_bind)\n    \n    print('Data process complete')\nelse:\n    # Preprocessed ECFP data by hengck23\n    # Reference: https://www.kaggle.com/competitions/leash-BELKA/discussion/492846#2748792\n    train_ecfp_path = '/kaggle/input/leash-bio-processed-dataset/train.ecfp4.packed.npz'\n    train_bind_path = '/kaggle/input/leash-bio-processed-dataset/train.bind.npz'\n    test_path = '/kaggle/input/leash-bio-processed-dataset/test.ecfp4.packed.npz'\n\n    train_ecfp = np.load(train_ecfp_path)\n    train_bind = np.load(train_bind_path)\n    test_ecfp = np.load(test_path)\n\n    print('train_ecfp files: ')\n    print(train_ecfp.files)\n    print('--------------------')\n    print('train_bind files: ')\n    print(train_bind.files)\n    print('--------------------')\n    print('test_ecfp files: ')\n    print(test_ecfp.files)","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:09:54.239293Z","iopub.execute_input":"2024-06-30T23:09:54.239659Z","iopub.status.idle":"2024-06-30T23:12:16.503557Z","shell.execute_reply.started":"2024-06-30T23:09:54.239628Z","shell.execute_reply":"2024-06-30T23:12:16.502557Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# For faster processing\ntest_path = '/kaggle/input/leash-bio-processed-dataset/test.ecfp4.packed.npz'\ntest_ecfp = np.load(test_path)","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:12:16.505356Z","iopub.execute_input":"2024-06-30T23:12:16.505771Z","iopub.status.idle":"2024-06-30T23:12:16.514563Z","shell.execute_reply.started":"2024-06-30T23:12:16.505731Z","shell.execute_reply":"2024-06-30T23:12:16.513791Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"if CFG.PREPROCESS:\n    train = np.array([np.concatenate((x, y)) for x,y in zip(train_ecfp, train_bind)])\nelse:\n    train = np.array([np.concatenate((x, y)) for x,y in zip(train_ecfp['ecfp'], train_bind['bind'])])\ntest = test_ecfp['ecfp']\n\ntrain.shape, test.shape","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:12:16.515588Z","iopub.execute_input":"2024-06-30T23:12:16.515823Z","iopub.status.idle":"2024-06-30T23:12:17.634623Z","shell.execute_reply.started":"2024-06-30T23:12:16.515799Z","shell.execute_reply":"2024-06-30T23:12:17.633716Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"bit = [f'bit{i}' for i in range(CFG.INP_LEN)]\ntrain = pd.DataFrame(train, columns=bit+['bind1', 'bind2', 'bind3'])\ntest = pd.DataFrame(test, columns=bit)","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:12:17.635906Z","iopub.execute_input":"2024-06-30T23:12:17.636325Z","iopub.status.idle":"2024-06-30T23:12:17.643031Z","shell.execute_reply.started":"2024-06-30T23:12:17.636290Z","shell.execute_reply":"2024-06-30T23:12:17.642111Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train.head()","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:12:17.644116Z","iopub.execute_input":"2024-06-30T23:12:17.644413Z","iopub.status.idle":"2024-06-30T23:12:17.667639Z","shell.execute_reply.started":"2024-06-30T23:12:17.644389Z","shell.execute_reply":"2024-06-30T23:12:17.666761Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test.head","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:12:17.668698Z","iopub.execute_input":"2024-06-30T23:12:17.668959Z","iopub.status.idle":"2024-06-30T23:12:17.838259Z","shell.execute_reply.started":"2024-06-30T23:12:17.668935Z","shell.execute_reply":"2024-06-30T23:12:17.837421Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Modeling: 1D Convolutional NN Architeture","metadata":{}},{"cell_type":"code","source":"def my_model():\n    with strategy.scope():\n        NUM_FILTERS = 32\n        hidden_dim = 128\n\n        inputs = tf.keras.layers.Input(shape=(CFG.INP_LEN,), dtype='int32')\n        x = tf.keras.layers.Embedding(input_dim=CFG.INP_LEN, output_dim=hidden_dim, input_length=CFG.INP_LEN, mask_zero = True)(inputs)\n        x = tf.keras.layers.Conv1D(filters=NUM_FILTERS, kernel_size=3,  activation='relu', padding='valid',  strides=1)(x)\n        x = tf.keras.layers.Conv1D(filters=NUM_FILTERS*2, kernel_size=3,  activation='relu', padding='valid',  strides=1)(x)\n        x = tf.keras.layers.Conv1D(filters=NUM_FILTERS*3, kernel_size=3,  activation='relu', padding='valid',  strides=1)(x)\n        x = tf.keras.layers.GlobalMaxPooling1D()(x)\n\n        x = tf.keras.layers.Dense(CFG.INP_LEN*4, activation='relu')(x)\n        x = tf.keras.layers.Dropout(0.1)(x)\n        x = tf.keras.layers.Dense(CFG.INP_LEN*4, activation='relu')(x)\n        x = tf.keras.layers.Dropout(0.1)(x)\n        x = tf.keras.layers.Dense(CFG.INP_LEN*2, activation='relu')(x)\n        x = tf.keras.layers.Dropout(0.1)(x)\n\n        outputs = tf.keras.layers.Dense(3, activation='sigmoid')(x)\n\n        model = tf.keras.models.Model(inputs=inputs, outputs=outputs)\n        optimizer = tf.keras.optimizers.AdamW(learning_rate=CFG.LR, weight_decay=CFG.WD)\n        loss = 'binary_crossentropy'\n        weighted_metrics = [tf.keras.metrics.AUC(curve='PR', name='avg_precision')]\n        model.compile(\n        loss=loss,\n        optimizer=optimizer,\n        weighted_metrics=weighted_metrics,\n        )\n        \n        return model\n    \nmodel = my_model()\n\nmodel.summary()\ntf.keras.backend.clear_session()","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:12:17.839349Z","iopub.execute_input":"2024-06-30T23:12:17.839628Z","iopub.status.idle":"2024-06-30T23:12:18.514599Z","shell.execute_reply.started":"2024-06-30T23:12:17.839601Z","shell.execute_reply":"2024-06-30T23:12:18.513651Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Data Training & Inference","metadata":{}},{"cell_type":"code","source":"from sklearn.model_selection import StratifiedKFold\nfrom sklearn.metrics import average_precision_score as APS\n\nFEATURES = test.columns\nTARGETS = ['bind1', 'bind2', 'bind3']\nskf = StratifiedKFold(n_splits=CFG.NBR_FOLDS, shuffle=True, random_state=42)\n\nall_preds = []\nfor fold,(train_idx, valid_idx) in enumerate(skf.split(train, train[TARGETS].sum(1))):\n    \n    if fold not in CFG.SELECTED_FOLDS:\n        continue;\n    \n    X_train = train.loc[train_idx, FEATURES]\n    y_train = train.loc[train_idx, TARGETS]\n    X_val = train.loc[valid_idx, FEATURES]\n    y_val = train.loc[valid_idx, TARGETS]\n\n    es = tf.keras.callbacks.EarlyStopping(patience=5, monitor=\"val_loss\", mode='min', verbose=1)\n    checkpoint = tf.keras.callbacks.ModelCheckpoint(monitor='val_loss', filepath=f\"model-{fold}.weights.h5\",\n                                                        save_best_only=True, save_weights_only=True, mode='min')\n    reduce_lr_loss = tf.keras.callbacks.ReduceLROnPlateau(monitor='val_loss', factor=0.05, patience=5, verbose=1)\n    model = my_model()\n    history = model.fit(\n            X_train, y_train,\n            validation_data=(X_val, y_val),\n            epochs=CFG.EPOCHS,\n            callbacks=[checkpoint, reduce_lr_loss, es],\n            batch_size=CFG.BATCH_SIZE,\n            verbose=1,\n        )\n    \n    model.load_weights(f\"model-{fold}.weights.h5\")\n    oof = model.predict(X_val, batch_size=CFG.BATCH_SIZE)\n    print('fold :', fold, 'CV score =', APS(y_val, oof, average = 'micro'))\n    \n    preds = model.predict(test, batch_size=CFG.BATCH_SIZE)\n    all_preds.append(preds)\n\npreds = np.mean(all_preds, 0)","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:12:18.518916Z","iopub.execute_input":"2024-06-30T23:12:18.519267Z","iopub.status.idle":"2024-06-30T23:12:50.179147Z","shell.execute_reply.started":"2024-06-30T23:12:18.519242Z","shell.execute_reply":"2024-06-30T23:12:50.178265Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"preds[:5]","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:12:50.180342Z","iopub.execute_input":"2024-06-30T23:12:50.180638Z","iopub.status.idle":"2024-06-30T23:12:50.187158Z","shell.execute_reply.started":"2024-06-30T23:12:50.180611Z","shell.execute_reply":"2024-06-30T23:12:50.186272Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Submission","metadata":{}},{"cell_type":"code","source":"tst = pd.read_parquet('/kaggle/input/leash-BELKA/test.parquet')\ntst_smile_unique = tst['molecule_smiles'].unique()\ntest_preds = np.column_stack((tst_smile_unique, preds))\n\npred_df = pd.DataFrame(test_preds, columns=['molecule_smiles', 'bind1', 'bind2', 'bind3'])\n\nresult = pd.merge(tst, pred_df, how='left', on='molecule_smiles')\n\nresult.head()","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:12:50.188502Z","iopub.execute_input":"2024-06-30T23:12:50.188779Z","iopub.status.idle":"2024-06-30T23:12:52.741169Z","shell.execute_reply.started":"2024-06-30T23:12:50.188752Z","shell.execute_reply":"2024-06-30T23:12:52.740107Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def get_bind(row):\n    if row['protein_name']=='BRD4':\n        return row['bind1']\n    elif row['protein_name']=='HSA':\n        return row['bind2']\n    else:\n        return row['bind3']\n\n\nsubmission = pd.DataFrame()\n\nsubmission['id'] = result['id']\nsubmission['binds'] = result.apply(get_bind, axis=1)\n\nsubmission.head()","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:12:52.742403Z","iopub.execute_input":"2024-06-30T23:12:52.742733Z","iopub.status.idle":"2024-06-30T23:13:22.411952Z","shell.execute_reply.started":"2024-06-30T23:12:52.742702Z","shell.execute_reply":"2024-06-30T23:13:22.410969Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"submission.to_csv('submission.csv', index = False)","metadata":{"execution":{"iopub.status.busy":"2024-06-30T23:13:22.413415Z","iopub.execute_input":"2024-06-30T23:13:22.414296Z","iopub.status.idle":"2024-06-30T23:13:27.460498Z","shell.execute_reply.started":"2024-06-30T23:13:22.414254Z","shell.execute_reply":"2024-06-30T23:13:27.459267Z"},"trusted":true},"execution_count":null,"outputs":[]}]}