{"cells":[{"metadata":{},"cell_type":"markdown","source":"## 加载功能"},{"metadata":{"trusted":true},"cell_type":"code","source":"%matplotlib inline\n\nimport pandas as pd\nimport numpy as np\n\nimport math\nimport gc\nimport copy\nimport os\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\n\nfrom lightgbm import LGBMRegressor","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#设置数据输入路径和导出路径\nDATA_PATH = '../input'\nSUBMISSIONS_PATH = './'\n#用原子序数代替原子名称\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":"#设置显示参数\npd.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":"#显示当前路径下的文件夹\nos.listdir(DATA_PATH)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#定义缩小消耗内存函数\ndef reduce_mem_usage(df, verbose=False):\n    numerics = ['int16', 'int32', 'int64', 'float16', 'float32', 'float64']\n    start_mem = df.memory_usage().sum() / 1024**2    \n    for col in df.columns:\n        col_type = df[col].dtypes\n        if col_type in numerics:\n            c_min = df[col].min()\n            c_max = df[col].max()\n            if str(col_type)[:3] == 'int':\n                if c_min > np.iinfo(np.int8).min and c_max < np.iinfo(np.int8).max:\n                    df[col] = df[col].astype(np.int8)\n                elif c_min > np.iinfo(np.int16).min and c_max < np.iinfo(np.int16).max:\n                    df[col] = df[col].astype(np.int16)\n                elif c_min > np.iinfo(np.int32).min and c_max < np.iinfo(np.int32).max:\n                    df[col] = df[col].astype(np.int32)\n                elif c_min > np.iinfo(np.int64).min and c_max < np.iinfo(np.int64).max:\n                    df[col] = df[col].astype(np.int64)  \n            else:\n                if c_min > np.finfo(np.float32).min and c_max < np.finfo(np.float32).max:\n                    df[col] = df[col].astype(np.float32)\n                else:\n                    df[col] = df[col].astype(np.float64)    \n    end_mem = df.memory_usage().sum() / 1024**2\n    if verbose: print('Mem. usage decreased to {:5.2f} Mb ({:.1f}% reduction)'.format(end_mem, 100 * (start_mem - end_mem) / start_mem))\n    return df","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 加载数据集"},{"metadata":{"trusted":true},"cell_type":"code","source":"train_csv = pd.read_csv('../input/data-of-distance-qm9-giba/fin_train.csv')\ntest_csv = pd.read_csv('../input/data-of-distance-qm9-giba/fin_test.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#减小消耗内存，防止内存不足\ntrain_csv = reduce_mem_usage(train_csv,verbose = True)\ntest_csv = reduce_mem_usage(test_csv,verbose = True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_csv.set_index('id',inplace = True)\ntest_csv.set_index('id',inplace = True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_csv.head(10)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_csv.head(10)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"structures_csv = pd.read_csv(f'{DATA_PATH}/champs-scalar-coupling/structures.csv')\n#使用新的数值型特征molecule_index代替字符型molecule_name\nstructures_csv['molecule_index'] = structures_csv.molecule_name.str.replace('dsgdb9nsd_', '').astype('int32')\n#删除molecule_name列\nstructures_csv = structures_csv[['molecule_index', 'atom_index', 'atom', 'x', 'y', 'z']]\n#使用原子序数代替原子名称\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 = reduce_mem_usage(structures_csv,verbose = True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#定义函数：使base和structures只保留coupling_type类型的数据\ndef 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":"#定义函数：根据molecule_index和atom_index_{index}合并base和structures\ndef 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":"#定义函数：根据molecule_index、atom_index_0、atom_index_1合并base和atoms\ndef 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":"#定义函数：根据molecule_index合并base和structures\ndef 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":"#定义函数：计算（x_0,y_0,z_0）和（x_1,y_1,z_1）两点之间的中点\ndef 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#定义函数：计算原子到中心的距离\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\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":"#定义函数：计算所有的点两两之间的距离\ndef 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":"#定义函数：计算structures中每个分子包含的原子数并将该数值合并到base中\ndef 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":"#定义函数：添加距离相关特征\ndef 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    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":"#定义函数：添加所有需要的特征标签\ndef take_n_atoms(df, n_atoms, four_start=4):\n    #除了距离相关的特征剩下的特征\n    labels = [\n        'linkM0',\n       'linkM1', 'inv_dist0R', 'mean_molecule_atom_0_dist_xyz',\n       'max_molecule_atom_0_dist_xyz', 'sd_molecule_atom_0_dist_xyz',\n       'mean_molecule_atom_1_dist_xyz', 'max_molecule_atom_1_dist_xyz',\n       'sd_molecule_atom_1_dist_xyz', 'coulomb_C.x', 'coulomb_H.x',\n       'coulomb_N.x', 'coulomb_O.x', 'yukawa_H.x', 'coulomb_C.y',\n       'coulomb_H.y', 'coulomb_O.y', 'yukawa_C.y', 'yukawa_H.y', 'vander_C.y',\n       'distC0', 'distC1', 'adH1', 'adC2', 'adC3', 'adC4',\n       'molecule_index', 'mulliken_min', 'mulliken_max', 'mulliken_atom_0',\n       'mulliken_atom_1'\n    ]\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":"#LGB模型参数\nLGB_PARAMS = {\n    'objective': 'regression',\n    'metric': 'mae',\n    'verbosity': -1,\n    'boosting_type': 'gbdt',\n    'learning_rate': 0.1455,\n    'num_leaves': 129,\n    'min_child_samples': 78,\n    'max_depth': 13,\n    'subsample_freq': 1,\n    'subsample': 0.88,\n    'bagging_seed': 15,\n    'reg_alpha': 0.10107001,\n    'reg_lambda': 0.300132,\n    'colsample_bytree': 1.0\n}","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 最终模型"},{"metadata":{"trusted":true},"cell_type":"code","source":"submission_csv = pd.read_csv(f'{DATA_PATH}/champs-scalar-coupling/sample_submission.csv', index_col='id')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#定义函数：创建训练数据和测试数据\ndef 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":"#定义函数：根据耦合的类型训练模型\ndef train_and_predict_for_one_coupling_type(coupling_type, submission, n_atoms, n_folds=5, 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 = LGBMRegressor(**LGB_PARAMS, n_estimators=8000, 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=1000, early_stopping_rounds=200)\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    plt.scatter(y_val, y_val_pred)   \n    plt.show()    \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}\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":"#每个耦合类型的CV得分\npd.DataFrame({'type': list(cv_scores.keys()), 'cv_score': list(cv_scores.values())})","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#CV得分均值\nnp.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.head(10)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"submission.to_csv(f'{SUBMISSIONS_PATH}/submission.csv')","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"codemirror_mode":{"name":"ipython","version":3},"file_extension":".py","mimetype":"text/x-python","name":"python","nbconvert_exporter":"python","pygments_lexer":"ipython3","version":"3.7.1"}},"nbformat":4,"nbformat_minor":1}