{"cells":[{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"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 in \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 \"../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# Any results you write to the current directory are saved as output.","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"I'm using https://www.kaggle.com/marcogorelli/criskiev-s-distances-more-estimators-groupkfold kernel as base, but added some features and finetuned lgb"},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","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 GroupKFold, train_test_split\nfrom sklearn.metrics import mean_absolute_error\n\nimport matplotlib.pyplot as plt\nimport seaborn as sns\n\n\nfrom lightgbm import LGBMRegressor\nimport lightgbm as lgb","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"DATA_PATH = '/kaggle/input/champs-scalar-coupling/'\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":"pd.set_option('display.max_colwidth', -1)\npd.set_option('display.max_rows', 120)\npd.set_option('display.max_columns', 120)","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":"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":"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\n#did find how to deal with data processing funtion at the moment of creation this, therefore another one center metween two atoms\ndef add_newcenter(df):\n    df['centx'] = ((df['x_1'] + df['x_0']) * np.float32(0.5))\n    df['centy'] = ((df['y_1'] + df['y_0']) * np.float32(0.5))\n    df['centz'] = ((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}'] +1e-10)**np.float32(2) +\n        (df[f'y_{suffix1}'] - df[f'y_{suffix2}'] +1e-10)**np.float32(2) + \n        (df[f'z_{suffix1}'] - df[f'z_{suffix2}'] +1e-10)**np.float32(2)\n    )**np.float32(0.5))\n\n# distance btw center of two target atoms and center of another non-target pair of atoms\ndef add_distance_ctoc(df, suffix1, suffix2):\n    df[f'dc_{suffix1}_{suffix2}'] = ((\n        ((df[f'x_{suffix1}'] + df[f'x_{suffix2}'])*np.float32(0.5) - df['centx'])**np.float32(2) +\n        ((df[f'y_{suffix1}'] + df[f'y_{suffix2}'])*np.float32(0.5) - df['centy'])**np.float32(2) + \n        ((df[f'z_{suffix1}'] + df[f'z_{suffix2}'])*np.float32(0.5) - df['centz'])**np.float32(2)\n    )**np.float32(0.5))\n    \n# cube distance is from some big book about nuclear or quantum physics\n# should try to add same for target atom vs non-target atom\ndef add_cube_dist(df, suffix1, suffix2):\n    df[f'd3_{suffix1}_{suffix2}'] = 1/((((\n        (df[f'x_{suffix1}'] - df[f'x_{suffix2}'] +1e-10)**np.float32(2) +\n        (df[f'y_{suffix1}'] - df[f'y_{suffix2}'] +1e-10)**np.float32(2) + \n        (df[f'z_{suffix1}'] - df[f'z_{suffix2}'] +1e-10)**np.float32(2)\n    )**np.float32(0.5)))**np.float32(3))","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)\n            add_distance_ctoc(df, i, vi)\n            add_cube_dist(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    \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    add_newcenter(base)    \n    add_center(atoms)\n\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\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\n    add_distances(full)\n\n    full['numofatoms'] = n_atoms\n    full['numofatoms'] = full['numofatoms'].astype('int8')\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            labels.append(f'd3_{i}_{j}')\n            labels.append(f'dc_{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\n# lgb.Dataset(data=)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"gc.collect()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"LGB_PARAMS = {\n    'objective': 'regression',\n    'metric': 'mae',\n    'verbosity': -1,\n    'boosting_type': 'gbdt',\n    'learning_rate': 0.07,\n    'num_leaves': 256,\n    'min_child_samples': 79,\n    'max_depth': 18,\n    'subsample_freq': 1,\n    'subsample': 0.9,\n    'bagging_seed': 11,\n    'reg_alpha': 0.1,\n    'reg_lambda': 0.3,\n#     'colsample_bytree': 1.0 #v.3.1\n    'colsample_bytree': 0.8 # v3.5\n}","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, full.molecule_index.values","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, full.molecule_index.values","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=5, n_splits=5, random_state=42):\n    print(f'*** Training Model for {coupling_type} ***')\n    \n    X_data, y_data, groups = 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    kfold = GroupKFold(n_splits=n_folds)\n\n    for fold, (train_index, val_index) in enumerate(kfold.split(X_data, y_data, groups=groups)):\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\n        model = LGBMRegressor(**LGB_PARAMS, n_estimators=15000, n_jobs = -1)\n        model.fit(X_train, y_train, \n            eval_set=[(X_train, y_train), (X_val, y_val)], eval_metric='mae',\n            verbose=500, early_stopping_rounds=300)\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        y_pred += model.predict(X_test)\n        \n        break\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":"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 = 7\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()))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"submission[submission['scalar_coupling_constant'] == 0].shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"submission.to_csv(f'{SUBMISSIONS_PATH}/chem_submission_lgbv35.csv')","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}