{"cells":[{"metadata":{},"cell_type":"markdown","source":"Credit to this kernel - Tried Catboost with Similar Distance method: https://www.kaggle.com/criskiev/distance-is-all-you-need-lb-1-481"},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"%matplotlib inline\n\nimport pandas as pd\nimport numpy as np\n\nimport math\nimport gc\nimport copy\n\nfrom sklearn.model_selection import KFold, train_test_split\nfrom sklearn.metrics import mean_absolute_error\n\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nfrom catboost import CatBoostRegressor\nfrom lightgbm import LGBMRegressor","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"DATA_PATH = '../input'\nSUBMISSIONS_PATH = './'\n# use atomic numbers to recode atomic names\nATOMIC_NUMBERS = {\n    'H': 1,\n    'C': 6,\n    'N': 7,\n    'O': 8,\n    'F': 9\n}","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_dtypes = {\n    'molecule_name': 'category',\n    'atom_index_0': 'int8',\n    'atom_index_1': 'int8',\n    'type': 'category',\n    'scalar_coupling_constant': 'float32'\n}\ntrain_csv = pd.read_csv(f'{DATA_PATH}/train.csv', index_col='id', dtype=train_dtypes)\ntrain_csv['molecule_index'] = train_csv.molecule_name.str.replace('dsgdb9nsd_', '').astype('int32')\ntrain_csv = train_csv[['molecule_index', 'atom_index_0', 'atom_index_1', 'type', 'scalar_coupling_constant']]\ntrain_csv.head(10)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print('Shape: ', train_csv.shape)\nprint('Total: ', train_csv.memory_usage().sum())\ntrain_csv.memory_usage()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_csv = pd.read_csv(f'{DATA_PATH}/test.csv', index_col='id', dtype=train_dtypes)\ntest_csv['molecule_index'] = test_csv['molecule_name'].str.replace('dsgdb9nsd_', '').astype('int32')\ntest_csv = test_csv[['molecule_index', 'atom_index_0', 'atom_index_1', 'type']]\ntest_csv.head(10)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"structures_dtypes = {\n    'molecule_name': 'category',\n    'atom_index': 'int8',\n    'atom': 'category',\n    'x': 'float32',\n    'y': 'float32',\n    'z': 'float32'\n}\nstructures_csv = pd.read_csv(f'{DATA_PATH}/structures.csv', dtype=structures_dtypes)\nstructures_csv['molecule_index'] = structures_csv.molecule_name.str.replace('dsgdb9nsd_', '').astype('int32')\nstructures_csv = structures_csv[['molecule_index', 'atom_index', 'atom', 'x', 'y', 'z']]\nstructures_csv['atom'] = structures_csv['atom'].replace(ATOMIC_NUMBERS).astype('int8')\nstructures_csv.head(10)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print('Shape: ', structures_csv.shape)\nprint('Total: ', structures_csv.memory_usage().sum())\nstructures_csv.memory_usage()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def build_type_dataframes(base, structures, coupling_type):\n    base = base[base['type'] == coupling_type].drop('type', axis=1).copy()\n    base = base.reset_index()\n    base['id'] = base['id'].astype('int32')\n    structures = structures[structures['molecule_index'].isin(base['molecule_index'])]\n    return base, structures","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def add_coordinates(base, structures, index):\n    df = pd.merge(base, structures, how='inner',\n                  left_on=['molecule_index', f'atom_index_{index}'],\n                  right_on=['molecule_index', 'atom_index']).drop(['atom_index'], axis=1)\n    df = df.rename(columns={\n        'atom': f'atom_{index}',\n        'x': f'x_{index}',\n        'y': f'y_{index}',\n        'z': f'z_{index}'\n    })\n    return df","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def add_atoms(base, atoms):\n    df = pd.merge(base, atoms, how='inner',\n                  on=['molecule_index', 'atom_index_0', 'atom_index_1'])\n    return df","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def merge_all_atoms(base, structures):\n    df = pd.merge(base, structures, how='left',\n                  left_on=['molecule_index'],\n                  right_on=['molecule_index'])\n    df = df[(df.atom_index_0 != df.atom_index) & (df.atom_index_1 != df.atom_index)]\n    return df","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def add_center(df):\n    df['x_c'] = ((df['x_1'] + df['x_0']) * np.float32(0.5))\n    df['y_c'] = ((df['y_1'] + df['y_0']) * np.float32(0.5))\n    df['z_c'] = ((df['z_1'] + df['z_0']) * np.float32(0.5))\n\ndef add_distance_to_center(df):\n    df['d_c'] = ((\n        (df['x_c'] - df['x'])**np.float32(2) +\n        (df['y_c'] - df['y'])**np.float32(2) + \n        (df['z_c'] - df['z'])**np.float32(2)\n    )**np.float32(0.5))\n\ndef add_distance_between(df, suffix1, suffix2):\n    df[f'd_{suffix1}_{suffix2}'] = ((\n        (df[f'x_{suffix1}'] - df[f'x_{suffix2}'])**np.float32(2) +\n        (df[f'y_{suffix1}'] - df[f'y_{suffix2}'])**np.float32(2) + \n        (df[f'z_{suffix1}'] - df[f'z_{suffix2}'])**np.float32(2)\n    )**np.float32(0.5))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def add_distances(df):\n    n_atoms = 1 + max([int(c.split('_')[1]) for c in df.columns if c.startswith('x_')])\n    \n    for i in range(1, n_atoms):\n        for vi in range(min(4, i)):\n            add_distance_between(df, i, vi)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def add_n_atoms(base, structures):\n    dfs = structures['molecule_index'].value_counts().rename('n_atoms').to_frame()\n    return pd.merge(base, dfs, left_on='molecule_index', right_index=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def build_couple_dataframe(some_csv, structures_csv, coupling_type, n_atoms=10):\n    base, structures = build_type_dataframes(some_csv, structures_csv, coupling_type)\n    base = add_coordinates(base, structures, 0)\n    base = add_coordinates(base, structures, 1)\n    \n    base = base.drop(['atom_0', 'atom_1'], axis=1)\n    atoms = base.drop('id', axis=1).copy()\n    if 'scalar_coupling_constant' in some_csv:\n        atoms = atoms.drop(['scalar_coupling_constant'], axis=1)\n        \n    add_center(atoms)\n    atoms = atoms.drop(['x_0', 'y_0', 'z_0', 'x_1', 'y_1', 'z_1'], axis=1)\n\n    atoms = merge_all_atoms(atoms, structures)\n    \n    add_distance_to_center(atoms)\n    \n    atoms = atoms.drop(['x_c', 'y_c', 'z_c', 'atom_index'], axis=1)\n    atoms.sort_values(['molecule_index', 'atom_index_0', 'atom_index_1', 'd_c'], inplace=True)\n    atom_groups = atoms.groupby(['molecule_index', 'atom_index_0', 'atom_index_1'])\n    atoms['num'] = atom_groups.cumcount() + 2\n    atoms = atoms.drop(['d_c'], axis=1)\n    atoms = atoms[atoms['num'] < n_atoms]\n\n    atoms = atoms.set_index(['molecule_index', 'atom_index_0', 'atom_index_1', 'num']).unstack()\n    atoms.columns = [f'{col[0]}_{col[1]}' for col in atoms.columns]\n    atoms = atoms.reset_index()\n    \n    # downcast back to int8\n    for col in atoms.columns:\n        if col.startswith('atom_'):\n            atoms[col] = atoms[col].fillna(0).astype('int8')\n            \n    atoms['molecule_index'] = atoms['molecule_index'].astype('int32')\n    \n    full = add_atoms(base, atoms)\n    add_distances(full)\n    \n    full.sort_values('id', inplace=True)\n    \n    return full","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def take_n_atoms(df, n_atoms, four_start=4):\n    labels = []\n    for i in range(2, n_atoms):\n        label = f'atom_{i}'\n        labels.append(label)\n\n    for i in range(n_atoms):\n        num = min(i, 4) if i < four_start else 4\n        for j in range(num):\n            labels.append(f'd_{i}_{j}')\n    if 'scalar_coupling_constant' in df:\n        labels.append('scalar_coupling_constant')\n    return df[labels]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"%%time\nfull = build_couple_dataframe(train_csv, structures_csv, '1JHN', n_atoms=10)\nprint(full.shape)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df = take_n_atoms(full, 7)\n# LightGBM performs better with 0-s then with NaN-s\ndf = df.fillna(0)\ndf.columns","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"X_data = df.drop(['scalar_coupling_constant'], axis=1).values.astype('float32')\ny_data = df['scalar_coupling_constant'].values.astype('float32')\n\nX_train, X_val, y_train, y_val = train_test_split(X_data, y_data, test_size=0.2, random_state=128)\nX_train.shape, X_val.shape, y_train.shape, y_val.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"model = CatBoostRegressor(iterations=750,\n                             learning_rate=0.2,\n                             depth=10,\n                             eval_metric='RMSE',\n                             random_seed = 42,\n                             bagging_temperature = 0.2,\n                             od_type='Iter',\n                             metric_period = 50,\n                             od_wait=20,\n                         reg_lambda=0.3)\nmodel.fit(X_train, y_train, \n        eval_set=[(X_train, y_train), (X_val, y_val)])\n\ny_pred = model.predict(X_val)\nnp.log(mean_absolute_error(y_val, y_pred))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cols = list(df.columns)\ncols.remove('scalar_coupling_constant')\ncols\ndf_importance = pd.DataFrame({'feature': cols, 'importance': model.feature_importances_})\nsns.barplot(x=\"importance\", y=\"feature\", data=df_importance.sort_values('importance', ascending=False));","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def build_x_y_data(some_csv, coupling_type, n_atoms):\n    full = build_couple_dataframe(some_csv, structures_csv, coupling_type, n_atoms=n_atoms)\n    \n    df = take_n_atoms(full, n_atoms)\n    df = df.fillna(0)\n    print(df.columns)\n    \n    if 'scalar_coupling_constant' in df:\n        X_data = df.drop(['scalar_coupling_constant'], axis=1).values.astype('float32')\n        y_data = df['scalar_coupling_constant'].values.astype('float32')\n    else:\n        X_data = df.values.astype('float32')\n        y_data = None\n    \n    return X_data, y_data","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def train_and_predict_for_one_coupling_type(coupling_type, submission, n_atoms, n_folds=3, n_splits=5, random_state=128):\n    print(f'*** Training Model for {coupling_type} ***')\n    \n    X_data, y_data = build_x_y_data(train_csv, coupling_type, n_atoms)\n    X_test, _ = build_x_y_data(test_csv, coupling_type, n_atoms)\n    y_pred = np.zeros(X_test.shape[0], dtype='float32')\n\n    cv_score = 0\n    \n    if n_folds > n_splits:\n        n_splits = n_folds\n    \n    kfold = KFold(n_splits=n_splits, shuffle=True, random_state=random_state)\n\n    for fold, (train_index, val_index) in enumerate(kfold.split(X_data, y_data)):\n        if fold >= n_folds:\n            break\n\n        X_train, X_val = X_data[train_index], X_data[val_index]\n        y_train, y_val = y_data[train_index], y_data[val_index]\n\n        model =  CatBoostRegressor(iterations=500,\n                             learning_rate=0.2,\n                             depth=10,\n                             eval_metric='RMSE',\n                             random_seed = 42,\n                             bagging_temperature = 0.2,\n                             od_type='Iter',\n                             metric_period = 50,\n                             od_wait=20,\n                         reg_lambda=0.3)\n        model.fit(X_train, y_train, \n            eval_set=[(X_train, y_train), (X_val, y_val)])\n\n        y_val_pred = model.predict(X_val)\n        val_score = np.log(mean_absolute_error(y_val, y_val_pred))\n        print(f'{coupling_type} Fold {fold}, logMAE: {val_score}')\n        \n        cv_score += val_score / n_folds\n        y_pred += model.predict(X_test) / n_folds\n        \n        \n    submission.loc[test_csv['type'] == coupling_type, 'scalar_coupling_constant'] = y_pred\n    return cv_score","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"submission_csv = pd.read_csv(f'{DATA_PATH}/sample_submission.csv', index_col='id')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"model_params = {\n    '1JHN': 7,\n    '1JHC': 10,\n    '2JHH': 9,\n    '2JHN': 9,\n    '2JHC': 9,\n    '3JHH': 9,\n    '3JHC': 10,\n    '3JHN': 10\n}\nN_FOLDS = 5\nsubmission = submission_csv.copy()\n\ncv_scores = {}\nfor coupling_type in model_params.keys():\n    cv_score = train_and_predict_for_one_coupling_type(\n        coupling_type, submission, n_atoms=model_params[coupling_type], n_folds=N_FOLDS)\n    cv_scores[coupling_type] = cv_score","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"pd.DataFrame({'type': list(cv_scores.keys()), 'cv_score': list(cv_scores.values())})","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"np.mean(list(cv_scores.values()))\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"submission[submission['scalar_coupling_constant'] == 0].shape\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"submission.head(10)\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"submission.to_csv(f'{SUBMISSIONS_PATH}/submission.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"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":4,"nbformat_minor":1}