{"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":50160,"databundleVersionId":7921029,"sourceType":"competition"}],"dockerImageVersionId":30665,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import sys\nfrom pathlib import Path\nimport subprocess\nimport os\nimport gc\nfrom glob import glob\n\nimport numpy as np\nimport pandas as pd\nimport polars as pl\nfrom datetime import datetime\nimport seaborn as sns\nimport matplotlib.pyplot as plt\n\nimport optuna\nfrom optuna.visualization import plot_param_importances\nfrom bayes_opt import BayesianOptimization\n\nimport warnings\nwarnings.filterwarnings('ignore')\n\nROOT = '/kaggle/input/home-credit-credit-risk-model-stability'","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2024-03-22T15:37:24.601956Z","iopub.execute_input":"2024-03-22T15:37:24.602368Z","iopub.status.idle":"2024-03-22T15:37:28.198577Z","shell.execute_reply.started":"2024-03-22T15:37:24.602327Z","shell.execute_reply":"2024-03-22T15:37:28.197316Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.model_selection import TimeSeriesSplit, GroupKFold, StratifiedGroupKFold\nfrom sklearn.base import BaseEstimator, RegressorMixin\nfrom sklearn.metrics import roc_auc_score\nimport lightgbm as lgb","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:37:28.200939Z","iopub.execute_input":"2024-03-22T15:37:28.201744Z","iopub.status.idle":"2024-03-22T15:37:29.112503Z","shell.execute_reply.started":"2024-03-22T15:37:28.201701Z","shell.execute_reply":"2024-03-22T15:37:29.111346Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class Pipeline:\n\n    def set_table_dtypes(df):\n        for col in df.columns:\n            if col in [\"case_id\", \"WEEK_NUM\", \"num_group1\", \"num_group2\"]:\n                df = df.with_columns(pl.col(col).cast(pl.Int64))\n            elif col in [\"date_decision\"]:\n                df = df.with_columns(pl.col(col).cast(pl.Date))\n            elif col[-1] in (\"P\", \"A\"):\n                df = df.with_columns(pl.col(col).cast(pl.Float64))\n            elif col[-1] in (\"M\",):\n                df = df.with_columns(pl.col(col).cast(pl.String))\n            elif col[-1] in (\"D\",):\n                df = df.with_columns(pl.col(col).cast(pl.Date))\n        return df\n\n    def handle_dates(df):\n        for col in df.columns:\n            if col[-1] in (\"D\",):\n                df = df.with_columns(pl.col(col) - pl.col(\"date_decision\"))  #!!?\n                df = df.with_columns(pl.col(col).dt.total_days()) # t - t-1\n        df = df.drop(\"date_decision\", \"MONTH\")\n        return df\n\n    def filter_cols(df):\n        for col in df.columns:\n            if col not in [\"target\", \"case_id\", \"WEEK_NUM\"]:\n                isnull = df[col].is_null().mean()\n                if isnull > 0.95:\n                    df = df.drop(col)\n        \n        for col in df.columns:\n            if (col not in [\"target\", \"case_id\", \"WEEK_NUM\"]) & (df[col].dtype == pl.String):\n                freq = df[col].n_unique()\n                if (freq == 1) | (freq > 200):\n                    df = df.drop(col)\n        \n        return df","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:37:29.113884Z","iopub.execute_input":"2024-03-22T15:37:29.114427Z","iopub.status.idle":"2024-03-22T15:37:29.129082Z","shell.execute_reply.started":"2024-03-22T15:37:29.114397Z","shell.execute_reply":"2024-03-22T15:37:29.127938Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class Aggregator:\n    \n    def num_expr(df):\n        cols = [col for col in df.columns if col[-1] in (\"P\", \"A\")]\n        expr_max = [pl.max(col).alias(f\"max_{col}\") for col in cols]\n        return expr_max\n    \n    def date_expr(df):\n        cols = [col for col in df.columns if col[-1] in (\"D\")]\n        expr_max = [pl.max(col).alias(f\"max_{col}\") for col in cols]\n        return expr_max\n    \n    def str_expr(df):\n        cols = [col for col in df.columns if col[-1] in (\"M\",)]\n        expr_max = [pl.max(col).alias(f\"max_{col}\") for col in cols]\n        return expr_max\n    \n    def other_expr(df):\n        cols = [col for col in df.columns if col[-1] in (\"T\", \"L\")]\n        expr_max = [pl.max(col).alias(f\"max_{col}\") for col in cols]\n        return expr_max\n    \n    def count_expr(df):\n        cols = [col for col in df.columns if \"num_group\" in col]\n        expr_max = [pl.max(col).alias(f\"max_{col}\") for col in cols]  # max & replace col name\n        return expr_max\n    \n    def get_exprs(df):\n        exprs = Aggregator.num_expr(df) + \\\n                Aggregator.date_expr(df) + \\\n                Aggregator.str_expr(df) + \\\n                Aggregator.other_expr(df) + \\\n                Aggregator.count_expr(df)\n\n        return exprs","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:37:29.132415Z","iopub.execute_input":"2024-03-22T15:37:29.133417Z","iopub.status.idle":"2024-03-22T15:37:29.1484Z","shell.execute_reply.started":"2024-03-22T15:37:29.13337Z","shell.execute_reply":"2024-03-22T15:37:29.147242Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def read_file(path, depth=None):\n    df = pl.read_parquet(path)\n    df = df.pipe(Pipeline.set_table_dtypes)\n    if depth in [1,2]:\n        df = df.group_by(\"case_id\").agg(Aggregator.get_exprs(df)) \n    return df\n\ndef read_files(regex_path, depth=None):\n    chunks = []\n    \n    for path in glob(str(regex_path)):\n        df = pl.read_parquet(path)\n        df = df.pipe(Pipeline.set_table_dtypes)\n        if depth in [1, 2]:\n            df = df.group_by(\"case_id\").agg(Aggregator.get_exprs(df))\n        chunks.append(df)\n    \n    df = pl.concat(chunks, how=\"vertical_relaxed\")\n    df = df.unique(subset=[\"case_id\"])\n    return df","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:37:29.149974Z","iopub.execute_input":"2024-03-22T15:37:29.150397Z","iopub.status.idle":"2024-03-22T15:37:29.164061Z","shell.execute_reply.started":"2024-03-22T15:37:29.150354Z","shell.execute_reply":"2024-03-22T15:37:29.163147Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def feature_eng(df_base, depth_0, depth_1, depth_2):\n    df_base = (\n        df_base\n        .with_columns(\n            month_decision = pl.col(\"date_decision\").dt.month(),\n            weekday_decision = pl.col(\"date_decision\").dt.weekday(),\n        )\n    )\n    for i, df in enumerate(depth_0 + depth_1 + depth_2):\n        df_base = df_base.join(df, how=\"left\", on=\"case_id\", suffix=f\"_{i}\")\n    df_base = df_base.pipe(Pipeline.handle_dates)\n    return df_base","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:37:29.165559Z","iopub.execute_input":"2024-03-22T15:37:29.166082Z","iopub.status.idle":"2024-03-22T15:37:29.177572Z","shell.execute_reply.started":"2024-03-22T15:37:29.166054Z","shell.execute_reply":"2024-03-22T15:37:29.17548Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def to_pandas(df_data, cat_cols=None):\n    df_data = df_data.to_pandas()\n    if cat_cols is None:\n        cat_cols = list(df_data.select_dtypes(\"object\").columns)\n    df_data[cat_cols] = df_data[cat_cols].astype(\"category\")\n    return df_data, cat_cols","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:37:29.17902Z","iopub.execute_input":"2024-03-22T15:37:29.179409Z","iopub.status.idle":"2024-03-22T15:37:29.191555Z","shell.execute_reply.started":"2024-03-22T15:37:29.179373Z","shell.execute_reply":"2024-03-22T15:37:29.190434Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ROOT            = Path(\"/kaggle/input/home-credit-credit-risk-model-stability\")\n\nTRAIN_DIR       = ROOT / \"parquet_files\" / \"train\"\nTEST_DIR        = ROOT / \"parquet_files\" / \"test\"","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:37:29.193205Z","iopub.execute_input":"2024-03-22T15:37:29.193587Z","iopub.status.idle":"2024-03-22T15:37:29.201301Z","shell.execute_reply.started":"2024-03-22T15:37:29.193546Z","shell.execute_reply":"2024-03-22T15:37:29.200202Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data_store = {\n    \"df_base\": read_file(TRAIN_DIR / \"train_base.parquet\"),\n    \"depth_0\": [\n        read_file(TRAIN_DIR / \"train_static_cb_0.parquet\"),\n        read_files(TRAIN_DIR / \"train_static_0_*.parquet\"),\n    ],\n    \"depth_1\": [\n        read_files(TRAIN_DIR / \"train_applprev_1_*.parquet\", 1),\n        read_file(TRAIN_DIR / \"train_tax_registry_a_1.parquet\", 1),\n        read_file(TRAIN_DIR / \"train_tax_registry_b_1.parquet\", 1),\n        read_file(TRAIN_DIR / \"train_tax_registry_c_1.parquet\", 1),\n        read_files(TRAIN_DIR / \"train_credit_bureau_a_1_*.parquet\", 1),\n        read_file(TRAIN_DIR / \"train_credit_bureau_b_1.parquet\", 1),\n        read_file(TRAIN_DIR / \"train_other_1.parquet\", 1),\n        read_file(TRAIN_DIR / \"train_person_1.parquet\", 1),\n        read_file(TRAIN_DIR / \"train_deposit_1.parquet\", 1),\n        read_file(TRAIN_DIR / \"train_debitcard_1.parquet\", 1),\n    ],\n    \"depth_2\": [\n        read_file(TRAIN_DIR / \"train_credit_bureau_b_2.parquet\", 2),\n        read_files(TRAIN_DIR / \"train_credit_bureau_a_2_*.parquet\", 2),\n    ]\n}","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:37:29.202772Z","iopub.execute_input":"2024-03-22T15:37:29.203271Z","iopub.status.idle":"2024-03-22T15:40:23.311392Z","shell.execute_reply.started":"2024-03-22T15:37:29.203194Z","shell.execute_reply":"2024-03-22T15:40:23.310139Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_train = feature_eng(**data_store)\nprint(\"train data shape:\\t\", df_train.shape)","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:40:23.317338Z","iopub.execute_input":"2024-03-22T15:40:23.318129Z","iopub.status.idle":"2024-03-22T15:40:34.801555Z","shell.execute_reply.started":"2024-03-22T15:40:23.318083Z","shell.execute_reply":"2024-03-22T15:40:34.800456Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data_store = {\n    \"df_base\": read_file(TEST_DIR / \"test_base.parquet\"),\n    \"depth_0\": [\n        read_file(TEST_DIR / \"test_static_cb_0.parquet\"),\n        read_files(TEST_DIR / \"test_static_0_*.parquet\"),\n    ],\n    \"depth_1\": [\n        read_files(TEST_DIR / \"test_applprev_1_*.parquet\", 1),\n        read_file(TEST_DIR / \"test_tax_registry_a_1.parquet\", 1),\n        read_file(TEST_DIR / \"test_tax_registry_b_1.parquet\", 1),\n        read_file(TEST_DIR / \"test_tax_registry_c_1.parquet\", 1),\n        read_files(TEST_DIR / \"test_credit_bureau_a_1_*.parquet\", 1),\n        read_file(TEST_DIR / \"test_credit_bureau_b_1.parquet\", 1),\n        read_file(TEST_DIR / \"test_other_1.parquet\", 1),\n        read_file(TEST_DIR / \"test_person_1.parquet\", 1),\n        read_file(TEST_DIR / \"test_deposit_1.parquet\", 1),\n        read_file(TEST_DIR / \"test_debitcard_1.parquet\", 1),\n    ],\n    \"depth_2\": [\n        read_file(TEST_DIR / \"test_credit_bureau_b_2.parquet\", 2),\n        read_files(TEST_DIR / \"test_credit_bureau_a_2_*.parquet\", 2),\n    ]\n}","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:40:34.802965Z","iopub.execute_input":"2024-03-22T15:40:34.803398Z","iopub.status.idle":"2024-03-22T15:40:35.637265Z","shell.execute_reply.started":"2024-03-22T15:40:34.80336Z","shell.execute_reply":"2024-03-22T15:40:35.636156Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_test = feature_eng(**data_store)\nprint(\"test data shape:\\t\", df_test.shape)","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:40:35.638681Z","iopub.execute_input":"2024-03-22T15:40:35.639047Z","iopub.status.idle":"2024-03-22T15:40:35.688104Z","shell.execute_reply.started":"2024-03-22T15:40:35.639017Z","shell.execute_reply":"2024-03-22T15:40:35.686921Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Drop the insignificant features\ndf_train = df_train.pipe(Pipeline.filter_cols)\ndf_test = df_test.select([col for col in df_train.columns if col != \"target\"])\n\nprint(\"train data shape:\\t\", df_train.shape)\nprint(\"test data shape:\\t\", df_test.shape)","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:40:35.689494Z","iopub.execute_input":"2024-03-22T15:40:35.689899Z","iopub.status.idle":"2024-03-22T15:40:39.107132Z","shell.execute_reply.started":"2024-03-22T15:40:35.689869Z","shell.execute_reply":"2024-03-22T15:40:39.105875Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_train, cat_cols = to_pandas(df_train)\ndf_test, cat_cols = to_pandas(df_test, cat_cols)","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:40:39.108573Z","iopub.execute_input":"2024-03-22T15:40:39.10896Z","iopub.status.idle":"2024-03-22T15:41:03.667138Z","shell.execute_reply.started":"2024-03-22T15:40:39.10893Z","shell.execute_reply":"2024-03-22T15:41:03.665922Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"del data_store\ngc.collect()","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:41:03.668612Z","iopub.execute_input":"2024-03-22T15:41:03.668984Z","iopub.status.idle":"2024-03-22T15:41:03.820131Z","shell.execute_reply.started":"2024-03-22T15:41:03.668954Z","shell.execute_reply":"2024-03-22T15:41:03.818691Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display(df_train.head())\ndisplay(df_test.head())","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:41:03.821601Z","iopub.execute_input":"2024-03-22T15:41:03.822637Z","iopub.status.idle":"2024-03-22T15:41:03.884353Z","shell.execute_reply.started":"2024-03-22T15:41:03.822602Z","shell.execute_reply":"2024-03-22T15:41:03.883181Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(\"Train is duplicated:\\t\", df_train[\"case_id\"].duplicated().any())\nprint(\"Train Week Range:\\t\", (df_train[\"WEEK_NUM\"].min(), df_train[\"WEEK_NUM\"].max()))\nprint()\nprint(\"Test is duplicated:\\t\", df_test[\"case_id\"].duplicated().any())\nprint(\"Test Week Range:\\t\", (df_test[\"WEEK_NUM\"].min(), df_test[\"WEEK_NUM\"].max()))\n\nsns.lineplot(\n    data=df_train,\n    x=\"WEEK_NUM\",\n    y=\"target\",\n)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:41:03.885686Z","iopub.execute_input":"2024-03-22T15:41:03.886126Z","iopub.status.idle":"2024-03-22T15:41:18.888042Z","shell.execute_reply.started":"2024-03-22T15:41:03.886086Z","shell.execute_reply":"2024-03-22T15:41:18.886892Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"base_train = df_train[[\"case_id\", \"WEEK_NUM\", \"target\"]].copy()","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:41:18.889436Z","iopub.execute_input":"2024-03-22T15:41:18.890251Z","iopub.status.idle":"2024-03-22T15:41:18.920809Z","shell.execute_reply.started":"2024-03-22T15:41:18.89021Z","shell.execute_reply":"2024-03-22T15:41:18.91963Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Extract the required columns as NumPy arrays\nweek_num = base_train[\"WEEK_NUM\"].values\ntarget = base_train[\"target\"].values\n# score = base_train[\"score\"].values","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:41:18.92215Z","iopub.execute_input":"2024-03-22T15:41:18.922566Z","iopub.status.idle":"2024-03-22T15:41:18.928522Z","shell.execute_reply.started":"2024-03-22T15:41:18.922529Z","shell.execute_reply":"2024-03-22T15:41:18.927383Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def gini_stability(score,w_m = 'w',complete=True, w_fallingrate=88.0, w_resstd=-0.5):\n    if w_m == 'w':\n#         gini_in_time = base.loc[:, [\"WEEK_NUM\", \"target\", \"score\"]]\\\n#                 .sort_values(\"WEEK_NUM\")\\\n#                 .groupby(\"WEEK_NUM\")[[\"target\", \"score\"]]\\\n#                 .apply(lambda x: 2*roc_auc_score(x[\"target\"], x[\"score\"])-1).tolist()\n\n\n        \n\n        # Sort the arrays based on week_num\n        sorted_indices = np.argsort(week_num)\n        week_num_sorted = week_num[sorted_indices]\n        target_sorted = target[sorted_indices]\n        score_sorted = score[sorted_indices]\n\n        # Calculate Gini in time using NumPy operations\n        gini_in_time = []\n        unique_weeks = np.unique(week_num_sorted)\n        for week in unique_weeks:\n            mask = week_num_sorted == week\n            auc_score = roc_auc_score(target_sorted[mask], score_sorted[mask])\n            gini_in_time.append(2 * auc_score - 1)\n\n        gini_in_time = np.array(gini_in_time)\n\n        x = np.arange(len(gini_in_time))\n        y = gini_in_time\n        a, b = np.polyfit(x, y, 1)\n        y_hat = a*x + b\n        residuals = y - y_hat\n        res_std = np.std(residuals)\n        avg_gini = np.mean(gini_in_time)\n        stability_score = avg_gini + w_fallingrate * min(0, a) + w_resstd * res_std\n        if complete:\n            plt.figure()\n            plt.scatter(x,y)\n            plt.plot(x,y_hat)\n            return stability_score,avg_gini,a,res_std\n        else:\n            return stability_score\n\n    elif w_m == 'm':\n        gini_in_time = base.loc[:, [\"MONTH\", \"target\", \"score\"]]\\\n                .sort_values(\"MONTH\")\\\n                .groupby(\"MONTH\")[[\"target\", \"score\"]]\\\n                .apply(lambda x: 2*roc_auc_score(x[\"target\"], x[\"score\"])-1).tolist()\n\n        x = np.arange(len(gini_in_time))\n        y = gini_in_time\n        a, b = np.polyfit(x, y, 1)\n        y_hat = a*x + b\n        residuals = y - y_hat\n        res_std = np.std(residuals)\n        avg_gini = np.mean(gini_in_time)\n        stability_score = avg_gini + w_fallingrate * min(0, a) + w_resstd * res_std\n        if complete:\n            plt.figure()\n            plt.scatter(x,y)\n            plt.plot(x,y_hat)\n            return stability_score,avg_gini,a,res_std\n        else:\n            return stability_score\n    else:\n        print(\"ERROR\")","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:41:18.930174Z","iopub.execute_input":"2024-03-22T15:41:18.931047Z","iopub.status.idle":"2024-03-22T15:41:18.947615Z","shell.execute_reply.started":"2024-03-22T15:41:18.931004Z","shell.execute_reply":"2024-03-22T15:41:18.946685Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X = df_train.drop(columns=[\"target\", \"case_id\", \"WEEK_NUM\"])\ny = df_train[\"target\"]\nweeks = df_train[\"WEEK_NUM\"]\n\ncv = StratifiedGroupKFold(n_splits=5, shuffle=False)","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:41:18.948964Z","iopub.execute_input":"2024-03-22T15:41:18.949883Z","iopub.status.idle":"2024-03-22T15:41:20.741827Z","shell.execute_reply.started":"2024-03-22T15:41:18.949826Z","shell.execute_reply":"2024-03-22T15:41:20.740449Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X.shape","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:41:20.743359Z","iopub.execute_input":"2024-03-22T15:41:20.743728Z","iopub.status.idle":"2024-03-22T15:41:20.751011Z","shell.execute_reply.started":"2024-03-22T15:41:20.743698Z","shell.execute_reply":"2024-03-22T15:41:20.749745Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# del df_train\n\n# gc.collect()","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:41:20.752304Z","iopub.execute_input":"2024-03-22T15:41:20.752633Z","iopub.status.idle":"2024-03-22T15:41:20.760519Z","shell.execute_reply.started":"2024-03-22T15:41:20.752595Z","shell.execute_reply":"2024-03-22T15:41:20.759601Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# def custom_loss(y_pred):\n#     base_train['score'] = y_pred\n#     stability_score = gini_stability(base_train,w_m='w',complete=False)\n#     return stability_score","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:30:07.684502Z","iopub.execute_input":"2024-03-22T15:30:07.685503Z","iopub.status.idle":"2024-03-22T15:30:07.694905Z","shell.execute_reply.started":"2024-03-22T15:30:07.685467Z","shell.execute_reply":"2024-03-22T15:30:07.693922Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# import torch\n# import torch.nn as nn\n# import torch.optim as optim\n# # from torch.autograd import grad, gradgrad","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:30:07.696206Z","iopub.execute_input":"2024-03-22T15:30:07.696769Z","iopub.status.idle":"2024-03-22T15:30:07.705947Z","shell.execute_reply.started":"2024-03-22T15:30:07.696738Z","shell.execute_reply":"2024-03-22T15:30:07.704761Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# # Define custom gradient and hessian calculation using PyTorch\n# def custom_grad_hess(y_true, y_pred):\n#     y_true_torch = torch.tensor(y_true, dtype=torch.float32)\n#     y_pred_torch = torch.tensor(y_pred, dtype=torch.float32, requires_grad=True)\n\n#     # Calculate custom loss using PyTorch\n#     loss_tensor = torch.tensor(custom_loss(y_pred), dtype=torch.float32)\n    \n#     # Ensure y_pred_torch requires gradient\n#     y_pred_torch.requires_grad_(True)\n#     print(y_pred_torch.shape)\n#     print(loss_tensor.shape)\n#     # Calculate gradient\n#     grad = torch.autograd.grad(loss_tensor, y_pred_torch)\n    \n#     # Calculate Hessian (second derivative)\n#     grad_val = grad.detach()  # Detach to avoid tracking history\n#     hess = torch.autograd.grad(grad_val.sum(), y_pred_torch, create_graph=True)[0]\n    \n#     return grad.numpy(), hess.numpy()","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:30:07.707419Z","iopub.execute_input":"2024-03-22T15:30:07.708065Z","iopub.status.idle":"2024-03-22T15:30:07.717192Z","shell.execute_reply.started":"2024-03-22T15:30:07.708035Z","shell.execute_reply":"2024-03-22T15:30:07.71602Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# class CustomLoss(nn.Module):\n#     def __init__(self, base, w_fallingrate=88.0, w_resstd=-0.5):\n#         super(CustomLoss, self).__init__()\n#         self.base = base\n#         self.w_fallingrate = w_fallingrate\n#         self.w_resstd = w_resstd\n\n#     def forward(self, input_tensor):\n#         # Calculate your custom loss function using input_tensor (model predictions)\n#         # You can access self.base, self.w_fallingrate, and self.w_resstd within this function\n#         # Example calculation:\n#         gini_in_time = ...  # Calculate Gini in time based on self.base and input_tensor\n#         # Compute the components of your custom loss function\n#         avg_gini = torch.mean(gini_in_time)\n#         a, b = ...  # Calculate 'a' and 'b' parameters from your Gini calculation\n#         residuals = gini_in_time - (a * torch.arange(len(gini_in_time)) + b)\n#         res_std = torch.std(residuals)\n#         loss = avg_gini + self.w_fallingrate * torch.min(torch.tensor(0.0), a) + self.w_resstd * res_std\n#         return loss","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:30:07.718184Z","iopub.execute_input":"2024-03-22T15:30:07.718508Z","iopub.status.idle":"2024-03-22T15:30:07.736147Z","shell.execute_reply.started":"2024-03-22T15:30:07.718482Z","shell.execute_reply":"2024-03-22T15:30:07.735132Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# # Define custom objective function for LightGBM using PyTorch gradient and hessian\n# def custom_objective(y_true, y_pred):\n#     grad, hess = custom_grad_hess(y_true, y_pred)\n#     return grad, hess","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:30:07.740795Z","iopub.execute_input":"2024-03-22T15:30:07.7414Z","iopub.status.idle":"2024-03-22T15:30:07.749973Z","shell.execute_reply.started":"2024-03-22T15:30:07.741364Z","shell.execute_reply":"2024-03-22T15:30:07.748446Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from scipy.optimize import approx_fprime\nimport tensorflow as tf","metadata":{"execution":{"iopub.status.busy":"2024-03-22T15:41:20.761459Z","iopub.execute_input":"2024-03-22T15:41:20.761814Z","iopub.status.idle":"2024-03-22T15:41:34.690162Z","shell.execute_reply.started":"2024-03-22T15:41:20.761785Z","shell.execute_reply":"2024-03-22T15:41:34.688953Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Custom accuracy loss function\ndef custom_loss(score):\n    score = np.array(score)\n    w_fallingrate=88.0\n    w_resstd=-0.5\n    \n# Sort the arrays based on week_num\n    sorted_indices = np.argsort(week_num)\n    week_num_sorted = week_num[sorted_indices]\n    target_sorted = target[sorted_indices]\n    score_sorted = score[sorted_indices]\n\n    # Calculate Gini in time using NumPy operations\n    gini_in_time = []\n    unique_weeks = np.unique(week_num_sorted)\n    for week in unique_weeks:\n        mask = week_num_sorted == week\n        auc_score = roc_auc_score(target_sorted[mask], score_sorted[mask])\n        gini_in_time.append(2 * auc_score - 1)\n\n    gini_in_time = np.array(gini_in_time)\n\n    x = np.arange(len(gini_in_time))\n    yy = gini_in_time\n    a, b = np.polyfit(x, yy, 1)\n    y_hat = a*x + b\n    residuals = yy - y_hat\n    res_std = np.std(residuals)\n    avg_gini = np.mean(gini_in_time)\n    stability_score = avg_gini + w_fallingrate * min(0, a) + w_resstd * res_std    \n    print(f'stability_score = {stability_score}')\n    return tf.constant(1-stability_score)\n\n# Compute the Hessian using numerical differentiation\ndef hessian(func, x):\n    print('hess')\n    n = len(x)\n    hess = np.zeros((n, n))\n    for i in range(n):\n        for j in range(n):\n            print(f'{i}\\n{j}')\n            hess[i, j] = approx_fprime(x, lambda x: approx_fprime(x, func, epsilon=1e-6, args=(), kwargs={'x': x})[j], epsilon=1e-6)[i]\n    return hess\n\n# Custom gradient and hessian calculation for accuracy loss\ndef custom_grad_hess(y_true, y_pred):\n#     y_true_torch = torch.tensor(y_true, dtype=torch.float32)\n#     y_pred_torch = torch.tensor(y_pred, dtype=torch.float32, requires_grad=True)\n\n#     # Calculate custom accuracy loss using PyTorch\n#     loss = custom_loss(y_pred_torch)\n    \n    # Compute gradient using automatic differentiation\n    x = tf.constant(y_pred)\n    with tf.GradientTape() as tape:\n        tape.watch(x)  # Watch the variable x\n        yg = custom_loss(x)\n    grad = tape.gradient(yg, x)\n    \n    # Calculate the Hessian of my_func with respect to x\n    hessian = tf.hessians(custom_loss, y_pred)[0]  # hessians returns a list, we take the first element\n\n    return grad, hess\n\n# Custom objective function for LightGBM using accuracy loss\ndef custom_objective(y_true, y_pred):\n    grad, hess = custom_grad_hess(y_true, y_pred)\n    return grad, hess","metadata":{"execution":{"iopub.status.busy":"2024-03-22T16:03:23.303945Z","iopub.execute_input":"2024-03-22T16:03:23.304423Z","iopub.status.idle":"2024-03-22T16:03:23.327141Z","shell.execute_reply.started":"2024-03-22T16:03:23.30439Z","shell.execute_reply":"2024-03-22T16:03:23.325903Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# import numpy as np\n# from scipy.optimize import approx_fprime\n\n# # Define your custom function\n# def custom_func(x):\n#     return np.sum(x**2)\n\n# # Define a point at which to evaluate the gradient and Hessian\n# x = np.array([1.0, 2.0, 3.0])\n\n# # Compute the gradient using approx_fprime from SciPy\n# grad = approx_fprime(x, custom_func, epsilon=1e-6)\n\n# # Compute the Hessian using numerical differentiation\n# def hessian(func, x):\n#     n = len(x)\n#     hess = np.zeros((n, n))\n#     for i in range(n):\n#         for j in range(n):\n#             hess[i, j] = approx_fprime(x, lambda x: approx_fprime(x, func, epsilon=1e-6, args=(), kwargs={'x': x})[j], epsilon=1e-6)[i]\n#     return hess\n\n# hess = hessian(custom_func, x)","metadata":{"execution":{"iopub.status.busy":"2024-03-22T16:03:23.809213Z","iopub.execute_input":"2024-03-22T16:03:23.80961Z","iopub.status.idle":"2024-03-22T16:03:23.815442Z","shell.execute_reply.started":"2024-03-22T16:03:23.809582Z","shell.execute_reply":"2024-03-22T16:03:23.814264Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Create and train the LGBMClassifier with custom loss\nmodel = lgb.LGBMClassifier(objective=custom_objective, metric=\"None\", random_state=42,verbose=0)\nmodel.fit(X, y, eval_set=[(X, y)])","metadata":{"execution":{"iopub.status.busy":"2024-03-22T16:03:24.351263Z","iopub.execute_input":"2024-03-22T16:03:24.351661Z","iopub.status.idle":"2024-03-22T16:04:27.614793Z","shell.execute_reply.started":"2024-03-22T16:03:24.351633Z","shell.execute_reply":"2024-03-22T16:04:27.613032Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import numpy as np\nfrom sklearn.datasets import make_classification\nfrom sklearn.model_selection import train_test_split\nfrom lightgbm import LGBMClassifier\nfrom scipy.special import expit\n\n\n\n# Custom log loss function\ndef custom_log_loss(y_true, y_pred):\n    base_train['score'] = y_pred\n    stability_score = gini_stability(base_train,w_m='w',complete=False)\n#     eps = 1e-15\n#     y_pred = np.clip(y_pred, eps, 1 - eps)  # Clip predictions to avoid log(0)\n#     loss = -np.mean(y_true * np.log(y_pred) + (1 - y_true) * np.log(1 - y_pred))\n    return stability_score\n\n# Custom objective function for LightGBM\ndef custom_objective(y_true, y_pred):\n    grad = expit(y_true) - expit(y_pred)\n    hess = np.maximum(expit(y_pred) * (1 - expit(y_pred)), 1e-6)\n    return grad, hess\n\n# Create and train the LGBMClassifier with custom loss\nmodel = LGBMClassifier(objective=custom_objective, metric=\"None\", random_state=42)\nmodel.fit(X_train, y_train, eval_set=[(X_test, y_test)], early_stopping_rounds=10, verbose=0)\n\n# Evaluate the model using custom log loss\ny_pred_proba = model.predict_proba(X_test)[:, 1]\ncustom_loss = custom_log_loss(y_test, y_pred_proba)\nprint(f\"Custom Log Loss: {custom_loss:.4f}\")","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"i = 0\nbest_params_list = []\nbest_value_list = []\nmodel_no_list = []\nfor idx_train, idx_valid in cv.split(X, y, groups=weeks):\n    if i == 0:\n        X_train, y_train = X.iloc[idx_train], y.iloc[idx_train]\n        X_valid, y_valid = X.iloc[idx_valid], y.iloc[idx_valid]\n\n        def objective(trial):\n            params = {\n    #             'boosting_type': 'gbdt',\n    #             \"objective\": \"binary\",\n    #             'metric': 'auc',\n                'lambda_l1': trial.suggest_loguniform('lambda_l1', 0.001, 0.1),\n                'lambda_l2': trial.suggest_loguniform('lambda_l2', 0.001, 0.1),\n                'num_leaves': trial.suggest_int('num_leaves', 2, 100),\n                'max_depth': trial.suggest_int('max_depth', 2, 15),\n                'learning_rate': trial.suggest_loguniform('learning_rate', 0.01, 0.3),\n                'feature_fraction': trial.suggest_uniform('feature_fraction', 0.5, 0.8),\n                'bagging_fraction': trial.suggest_uniform('bagging_fraction', 0.5, 0.8),\n                'bagging_freq': trial.suggest_int('bagging_freq', 4, 10),\n                'min_child_samples': trial.suggest_int('min_child_samples', 5, 50),\n                'subsample': trial.suggest_uniform('subsample', 0.1, 1.0),\n                'colsample_bytree': trial.suggest_uniform('colsample_bytree', 0.5, 1.0),\n                'colsample_bynode': trial.suggest_uniform('colsample_bynode', 0.5, 1.0),\n                'min_split_gain': trial.suggest_loguniform('min_split_gain', 0.1, 0.5),\n#                 'min_data_in_bin': trial.suggest_int('min_data_in_bin', 5, 30),\n#                 'max_bin': trial.suggest_int('max_bin', 100, 500),\n    #             'extra_trees':True,\n    #             'n_estimators': 1000,\n    #             'device': 'gpu'\n    #             'random_state': 42,\n    #             'verbose': -1\n            }\n\n            model = lgb.LGBMClassifier(**params,boosting_type='gbdt',objective='binary',metric='auc',extra_trees=True,n_estimators=1000,device='gpu',random_state=42,verbose=-1)\n            model.fit(X_train, y_train)\n            y_pred_proba = model.predict_proba(X)[:, 1]\n        #         auc_score = roc_auc_score(y_train, y_pred_proba)\n            base_train[\"score\"] = y_pred_proba\n            stability_score = gini_stability(base_train,w_m='w',complete=False)\n            return stability_score\n\n        # Create a study object and optimize the objective function using Bayesian Optimization\n        study = optuna.create_study(direction='maximize')\n        study.optimize(objective, n_trials=25)\n\n        # Print the best hyperparameters and their corresponding AUC score\n        best_params = study.best_params\n        best_stability_score = study.best_value\n        print(\"Best Hyperparameters:\", best_params)\n        print(\"Best stability Score:\", best_stability_score)\n\n        model = lgb.LGBMClassifier(**best_params,boosting_type='gbdt',objective='binary',metric='auc',extra_trees=True,n_estimators=1000,device='gpu',random_state=42,verbose=-1)\n        model.fit(X_train, y_train)\n\n        model.booster_.save_model(f'lgbmmodel_{i}.txt')\n\n        best_params_list.append(study.best_params)\n        best_value_list.append(study.best_value)\n        model_no_list.append(i)\n    \n    i+=1\nmodel_best_params_df = pd.DataFrame(best_params_list)\nmodel_best_params_df['stability_value'] = best_value_list\nmodel_best_params_df['model no'] = model_no_list\nmodel_best_params_df.to_csv('modelBestParams.csv')","metadata":{"execution":{"iopub.status.busy":"2024-03-22T02:49:25.235023Z","iopub.execute_input":"2024-03-22T02:49:25.235844Z","iopub.status.idle":"2024-03-22T02:50:33.464478Z","shell.execute_reply.started":"2024-03-22T02:49:25.235805Z","shell.execute_reply":"2024-03-22T02:50:33.463317Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# X = df_train.drop(columns=[\"target\", \"case_id\", \"WEEK_NUM\"])\n# y = df_train[\"target\"]\n# weeks = df_train[\"WEEK_NUM\"]\n\n# cv = StratifiedGroupKFold(n_splits=5, shuffle=False,random_state=42)\n\n# params = {\n#     \"boosting_type\": \"gbdt\",\n#     \"objective\": \"binary\",\n#     \"metric\": \"auc\",\n#     \"max_depth\": 10,\n#     \"learning_rate\": 0.05,\n#     \"n_estimators\": 1000,\n#     \"colsample_bytree\": 0.8,\n#     \"colsample_bynode\": 0.8,\n#     \"verbose\": -1,\n#     \"random_state\": 42,\n#     \"reg_alpha\": 0.1,\n#     \"reg_lambda\": 10,\n#     \"extra_trees\":True,\n#     'num_leaves':64,\n#     \"device\": \"gpu\", \n#     \"verbose\": -1\n# }\n\n# fitted_models = []\n# cv_scores = []\n\n# for idx_train, idx_valid in cv.split(X, y, groups=weeks):\n#     X_train, y_train = X.iloc[idx_train], y.iloc[idx_train]\n#     X_valid, y_valid = X.iloc[idx_valid], y.iloc[idx_valid]\n    \n#     model = lgb.LGBMClassifier(**params)\n#     model.fit(\n#         X_train, y_train,\n#         eval_set = [(X_valid, y_valid)],\n#         callbacks = [lgb.log_evaluation(200), lgb.early_stopping(50)] )\n#     fitted_models.append(model)\n    \n#     y_pred_valid = model.predict_proba(X_valid)[:,1]\n#     auc_score = roc_auc_score(y_valid, y_pred_valid)\n#     cv_scores.append(auc_score)\n    \n# print(\"CV AUC scores: \", cv_scores)\n# print(\"Maximum CV AUC score: \", max(cv_scores))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# del df_train\n# del X\n# del y\n# del weeks\n# del X_train\n# del y_train\n# del X_valid\n# del y_valid\n\n# gc.collect()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Submission","metadata":{}},{"cell_type":"code","source":"# X_test = df_test.drop(columns=[\"WEEK_NUM\"])\n# X_test = X_test.set_index(\"case_id\")\n\n# # lgb_pred = pd.Series(model.predict_proba(X_test)[:, 1], index=X_test.index)\n# lgb_pred_l = [item.predict_proba(X_test)[:, 1] for item in fitted_models]\n# lgb_pred = np.mean(np.stack(lgb_pred_l),axis=0)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# df_subm = pd.read_csv(ROOT / \"sample_submission.csv\")\n# df_subm = df_subm.set_index(\"case_id\")\n\n# df_subm[\"score\"] = lgb_pred","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# df_subm.head()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# df_subm.to_csv(\"submission.csv\")","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# mst work this time","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}