{"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.10.14","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":81933,"databundleVersionId":9643020,"sourceType":"competition"}],"dockerImageVersionId":30761,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":false},"papermill":{"default_parameters":{},"duration":151.830354,"end_time":"2024-09-25T20:44:41.405417","environment_variables":{},"exception":null,"input_path":"__notebook__.ipynb","output_path":"__notebook__.ipynb","parameters":{},"start_time":"2024-09-25T20:42:09.575063","version":"2.6.0"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"**<font size=\"5\">[Child Mind Institute — Problematic Internet Use](https://www.kaggle.com/competitions/child-mind-institute-problematic-internet-use)</font>**\n\nAfter generating a large number of features from the actigraphy data (over 900), this notebook trains 2 Gradient Boosting Machines (GBMs) using XGBoost and LightGBM (CatBoost didn't work for me). Finally, it combines many versions of predictions of 1 model using a hard voting ensemble.\n\nThe dataset is highly noisy, which results in a high standard deviation in model performance metrics. To secure a strong position on the private leaderboard, it was essential to develop a robust pipeline that minimizes dependence on randomness. Specifically, I aimed for models with low standard deviation in theirs performances scores when varying the random_state parameter. Such models will have a stable QWK score on the public leaderboard when I change the random seed (see my comment in this [discussion](https://www.kaggle.com/competitions/child-mind-institute-problematic-internet-use/discussion/550074))\n\nTo reduce the high standard deviation of the score :\n* I applied stratification based on the target variable and the most important features.\n* I excluded actigraphy data from individuals who wore their accelerometers for less than a few days.\n* I removed outliers during the training phase (but kept them for validation).\n* I eliminated many irrelevant features, ensuring that no model used more than 15 features.\n* I combined multiple predictions using a hard voting ensemble.\n* and I submitted my final submission on the public leaderboard several times with different random states to ensure my score had the lowest possible variance.\n\n\nI choose GBMs because the dataset contains many missing (NaN) values, and I did not fully trust any imputation methods for handling them.  \nThe GBMs differ in several ways : they are trained on various sets of features, and their growth policies vary depending on the specific techniques employed by the three GBM frameworks.  \n\nActigraphy features was selected as follows :  \n* Forward Selection : I started by training a strong LGBM model without any actigraphy features and initialized an empty list of retained actigraphy features. Each actigraphy feature that improved the previous model’s performance was added to this list.\n* Backward Selection : Starting from a GBM model with all the actigraphy features from the previous list, I performed backward selection by iteratively removing the least useful features, one at a time. The least useful features were identified using the feature importance permutation technique.\n\nAnd I used optuna to fit optimal hyperparameters :\n* hyperparameters of GBM (max_depth, num_leaves, etc.)\n* number of outliers to drop,\n* minimum wearing number of days of accelerometer to consider actigraphy data could be usefull.\n\nCompared to the best notebook on the public leaderboard, I did not impute missing values, I did not use a neural network (to much missing values), I did not impute missing target values (!!!), I only kept records with a target for training, I did not blindly consider the 60 features from the autoencoder, and I did not use the autoencoder.\n\nI didn't use samples with missing target because I didn't know to use them : they had no actigraphy data !\n\nIn this way, I was looking for strong results, independent of randomness and luck.\n\n**Credits :**\n* [@Sheikh Muhammad Abdullah](https://www.kaggle.com/abdmental01) for parallel process_file function in https://www.kaggle.com/code/abdmental01/cmi-best-single-model\n* [@AmbrosM](https://www.kaggle.com/ambrosm) for regression idea instead of classification, for thresholds optimizer, for left_handed/right_handed features idea and for light features idea, in https://www.kaggle.com/code/ambrosm/piu-eda-which-makes-sense\n* [@M abdullah](https://www.kaggle.com/abdullah0a) for features ideas in https://www.kaggle.com/code/abdullah0a/ensamble-models\n* [@Daniel Dewey](https://www.kaggle.com/dan3dewey) for bed times features ad enmoXlight feature in https://www.kaggle.com/competitions/child-mind-institute-problematic-internet-use/discussion/551202\n\nVersions :   \n\n|Version|seed|Final CV|Public|Comment|  \n|--:|--:|:---:|:---:|---:|  \n| 53|1013|.4955|.447|XGB23 is not ready - 30 repeats|  \n| 54|3033|.4950|.441| - 30 repeats 27 * LGBM24 + 24 * xgb23|  \n| 55|4031|.4943|.438| - 30 repeats  23 * LGBM24 + 14 * xgb23|\n| 56|4037|.4888|.435| - 30 repeats 27 * lgbm24 + 16 lgbm23b|\n| 57|4137|.4949|.44| - 30 repeats 23 * LGBM23b + 10 * xgb23|\n| 58|4537|.4937|.436| - 30 repeats 27 * LGBM24 + 2 * xgb23|\n| 92|9090|.4928 ±.0009|.427| hard vote of 97 xgb24, xgb24>lgbm23c |\n| 93|6969|.4955 ±.0007|.434| hard vote of 89 xgb24, xgb24>lgbm23c |\n| 95|  86|.4951 ±.0018|.   | hard vote of 69 xgb24, xgb24>lgbm24 & 26 |\n| 98| 105|.4955 ±.0023|.   | hard vote of 51 lgbm24 > xgb24 & lgbm26 |\n| 99|3085|.4965 ±.0014|.   | hard vote of 69 xgb24, xgb24>lgbm24 & 26 |\n|100|3185|.4960 ±.0007|.   | hard vote of 91 lgbm24 > xgb24 & lgbm26 |\n|101|3495|.4957 ±.0007|.   | hard vote of 97 xgb24, xgb24>lgbm24 & 26 |\n|102|9425|.4927 ±.0007|.   | hard vote of 97 lgbm24 > xgb24 & lgbm26 |\n|103|3495|.4957 ±.0007|.434| hard vote of 97 xgb24, xgb24>lgbm24 & 26 |\n|104|   7|.4970 ±.0009|.442| hard vote of 97 lgbm24, lgbm24>xgb24 & 26 |\n|105|2025|.4972 ±.0006|.431| ==> hard vote of 97 xgb24>lgbm24 <==|\n|106/107|8048|.4956 ±.0023|.441| hard vote of 57 lgbm24, lgbm24>xgb27 & 24 |\n|108|3036|.4950 ±.0009|.| hard vote of 97 xgb24, xgb24>xgb27&lgbm24 |\n|109|2021|.4964 ±.0007|.| hard vote of 97 xgb24, xgb24>xgb27&lgbm24 |\n|110|21|.4967 ±.0011|.| hard vote of 97 xgb24, xgb24>xgb27&lgbm24 |\n|111|90|.4916 ±.0007|.| hard vote of 97 xgb24, xgb24>xgb27&lgbm24 |\n|112|2025|.4972 ±.0006|.431| ==> hard vote of 97 xgb24>lgbm24 <==|\n\n# Libraries & parameters & data","metadata":{"papermill":{"duration":0.009626,"end_time":"2024-09-25T20:42:13.099198","exception":false,"start_time":"2024-09-25T20:42:13.089572","status":"completed"},"tags":[]}},{"cell_type":"code","source":"n_repeats, n_splits, seed, keep_files, do_feat_imp = 10, 3, 1013, True, True \nn_repeats, n_splits, seed, keep_files, do_feat_imp = 150, 5, 1321, False, True \n#n_repeats, n_splits, seed, keep_files, do_feat_imp = 10, 5, 97, False, True \n\ndebug = False\nif n_splits == 3: debug = True","metadata":{"_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","papermill":{"duration":9.266873,"end_time":"2024-09-25T20:42:22.37517","exception":false,"start_time":"2024-09-25T20:42:13.108297","status":"completed"},"tags":[],"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T21:27:40.878458Z","iopub.execute_input":"2024-12-17T21:27:40.879013Z","iopub.status.idle":"2024-12-17T21:27:40.888176Z","shell.execute_reply.started":"2024-12-17T21:27:40.878964Z","shell.execute_reply":"2024-12-17T21:27:40.886918Z"},"_kg_hide-input":false},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import os\npath = \"/kaggle/input\" if os.path.isdir(\"/kaggle/input/\") else \".\"\n\ncuda = False\nif path == \"/kaggle/input\":\n    import torch\n    cuda = torch.cuda.is_available()\nprint(f\"Is GPU/CUDA available : {cuda}\")\n\noutput_path = \"output\"\nif not os.path.exists(output_path):\n    os.makedirs(output_path)\n\nimport numpy as np\nimport pandas as pd\npd.set_option('display.max_columns', 100)\npd.set_option('display.max_rows', 500)\n\nimport matplotlib.pyplot as plt\nimport seaborn as sns\n\nfrom concurrent.futures import ThreadPoolExecutor\nfrom tqdm import tqdm\nfrom copy import deepcopy\nimport warnings\n\nimport scipy.stats as stats\nfrom scipy.optimize import minimize, curve_fit\nfrom math import perm, comb, factorial\n\nfrom sklearn import set_config\nset_config(transform_output=\"pandas\")\n\nfrom sklearn.base import TransformerMixin, BaseEstimator\nfrom sklearn.pipeline import make_pipeline\nfrom sklearn.compose import ColumnTransformer\nfrom sklearn.preprocessing import OrdinalEncoder\n\nfrom sklearn.model_selection import RepeatedStratifiedKFold, StratifiedKFold\n\nfrom sklearn.impute import SimpleImputer\nfrom sklearn.ensemble import IsolationForest\n\nfrom sklearn.metrics import cohen_kappa_score\ndef score_(y_true, y_pred):\n    return cohen_kappa_score(y_true, y_pred, weights = 'quadratic')\n    \nfrom xgboost import XGBRegressor\nfrom lightgbm import LGBMRegressor, early_stopping, log_evaluation\nfrom catboost import CatBoostRegressor\n\nimport psutil\n__n_cores = psutil.cpu_count()     # Available CPU cores\nprint(f\"N CPU Cores : {__n_cores}\")\nfrom multiprocessing import Pool   # Multiprocess Runs\ndef df_parallelize_run(func, list_params):\n    num_cores = np.min([__n_cores, len(list_params)])\n    pool = Pool(num_cores)\n    res = pool.map(func, list_params)\n    pool.close()\n    pool.join()    \n    return res\n\nfrom functools import partial","metadata":{"_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","papermill":{"duration":9.266873,"end_time":"2024-09-25T20:42:22.37517","exception":false,"start_time":"2024-09-25T20:42:13.108297","status":"completed"},"tags":[],"trusted":true,"execution":{"iopub.status.busy":"2024-12-18T21:20:00.899262Z","iopub.execute_input":"2024-12-18T21:20:00.899677Z","iopub.status.idle":"2024-12-18T21:20:07.775408Z","shell.execute_reply.started":"2024-12-18T21:20:00.899639Z","shell.execute_reply":"2024-12-18T21:20:07.774326Z"},"_kg_hide-input":true,"jupyter":{"source_hidden":true}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"end = \"\\033[0m\" # reset\n\nbold       = \"\\033[1m\" ; resetbold       = \"\\033[21m\"\nunderline = \"\\033[4m\"  ; resetunderline  = \"\\033[24m\"\nblink      = \"\\033[5m\" ; resetblink      = \"\\033[25m\"\nreverse    = \"\\033[7m\" ; resetreverse    = \"\\033[27m\"\n\nDefault      = \"\\033[39m\" ; Black        = \"\\033[30m\" ; White        = \"\\033[97m\"\nRed          = \"\\033[31m\" ; LightRed     = \"\\033[91m\"\nGreen        = \"\\033[32m\" ; LightGreen   = \"\\033[92m\"\nYellow       = \"\\033[33m\" ; LightYellow  = \"\\033[93m\"\nBlue         = \"\\033[34m\" ; LightBlue    = \"\\033[94m\"\nMagenta      = \"\\033[35m\" ; LightMagenta = \"\\033[95m\"\nCyan         = \"\\033[36m\" ; LightCyan    = \"\\033[96m\"\nDarkGray     = \"\\033[90m\" ; LightGray    = \"\\033[37m\"\n\nBackgroundDefault     = \"\\033[49m\"  ; BackgroundBlack        = \"\\033[40m\"\nBackgroundRed         = \"\\033[41m\"  ; BackgroundLightRed     = \"\\033[101m\"\nBackgroundGreen       = \"\\033[42m\"  ; BackgroundLightGreen   = \"\\033[102m\"\nBackgroundYellow      = \"\\033[43m\"  ; BackgroundLightYellow  = \"\\033[103m\"\nBackgroundBlue        = \"\\033[44m\"  ; BackgroundLightBlue    = \"\\033[104m\"\nBackgroundMagenta     = \"\\033[45m\"  ; BackgroundLightMagenta = \"\\033[105m\"\nBackgroundCyan        = \"\\033[46m\"  ; BackgroundLightCyan    = \"\\033[106m\"\nBackgroundDarkGray    = \"\\033[100m\" ; BackgroundLightGray    = \"\\033[47m\"\n\nbold_blue = bold + LightBlue","metadata":{"_kg_hide-input":true,"papermill":{"duration":0.024825,"end_time":"2024-09-25T20:42:22.411123","exception":false,"start_time":"2024-09-25T20:42:22.386298","status":"completed"},"tags":[],"jupyter":{"source_hidden":true},"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T21:27:43.565299Z","iopub.execute_input":"2024-12-17T21:27:43.565719Z","iopub.status.idle":"2024-12-17T21:27:43.573961Z","shell.execute_reply.started":"2024-12-17T21:27:43.565683Z","shell.execute_reply":"2024-12-17T21:27:43.572830Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Model parameters","metadata":{}},{"cell_type":"code","source":"all_params = {\n    \n    # My best model\n    \"xgb24\":{ # XGB24 after some cleaning & Optuna V63 n°8\n        \"regressor\" : XGBRegressor,\n        \"params\"    : {\n            'learning_rate'         : .01,\n            'alpha'                 : 4.199995292665447, \n            'colsample_bytree'      : 0.8355226517382207, \n            'lambda'                : 18.73334313290196, \n            'max_bin'               : 185, \n            'max_depth'             : 5, \n            'min_child_weight'      : 26, \n            'subsample'             : 0.413527823331391,\n        },\n        \"features\"  : [\n            'Basic_Demos-Age', 'Basic_Demos-Sex', 'Physical-Height', 'FGC-FGC_PU', 'SDS-SDS_Total_Raw', 'PreInt_EduHx-computerinternet_hoursday', \n            'DEE_Weight', 'q95_max_light_school', 'norm_yz_evening_q90', 'anglez_evening_kurt', 'q95_med_enmo_evening', 'q75_med_norm_xy_all', \n        ],\n        'isolation'     : 0.029386967502224397, \n        \"ndays_min\"     : 2,\n        'target_mapper' : {0:0, 1:1, 2:2, 3:3},\n    },\n\n    # My second best model\n    \"lgbm24\" : { # Optuna v60 n°5\n        \"regressor\": LGBMRegressor,\n        \"params\"   : {\n            \"learning_rate\"         : .02,\n            'colsample_bytree': 0.679692229840442, \n            'min_child_samples': 139, \n            'num_leaves': 135, \n            'reg_alpha': 0.3599938818790582, \n            'reg_lambda': 13.010442911581617, \n            'subsample': 0.3567277990838326\n        },\n        \"features\" : [\n            'Basic_Demos-Age', 'Basic_Demos-Sex', 'Physical-Height', 'FGC-FGC_PU', 'SDS-SDS_Total_Raw', 'PreInt_EduHx-computerinternet_hoursday', \n            'q95_max_light_school', 'q90_med_norm_xy_afternoon', 'anglez_evening_kurt', 'norm_xz_evening_std', 'q90_max_light_wevening', \n            'norm_yz_night_q75', 'LST_TBW', \n        ],\n        \"isolation\"     : 0.04392243763264256, \n        \"ndays_min\"     : 2,\n        'target_mapper' : {0:0, 1:1, 2:2, 3:3},\n    },\n}\n\n\nall_params = {k:v for k,v in list(all_params.items()) if k[:1] != \"_\"}\nprint(f\"Number of models : {len(all_params)}\")\n\nfeatures_to_create  = []\nfor _, p in all_params.items():\n    features_to_create.extend([f for f in p[\"features\"] if f not in features_to_create])\nprint(features_to_create)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T21:27:54.137386Z","iopub.execute_input":"2024-12-17T21:27:54.137845Z","iopub.status.idle":"2024-12-17T21:27:54.156311Z","shell.execute_reply.started":"2024-12-17T21:27:54.137809Z","shell.execute_reply":"2024-12-17T21:27:54.155046Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## Actigraphy data & FE\nCredits : \n* [Sheikh Muhammad Abdullah](https://www.kaggle.com/abdmental01) for parallel process_file function in https://www.kaggle.com/code/abdmental01/cmi-best-single-model\n* [@Daniel Dewey](https://www.kaggle.com/dan3dewey) for bed times features ad enmoXlight feature in https://www.kaggle.com/competitions/child-mind-institute-problematic-internet-use/discussion/551202. But I didn't retain those features...\n\nI wasn't very original. I considered multiple instances of combinations between acceleration or light variables, with the time of day and the day of the week. It's a kind of brute force feature engineering.","metadata":{}},{"cell_type":"code","source":"features_ts_to_create = []\n\nfor f in [\"enmo\", \"X\", \"Y\", \"Z\", \"anglez\"]:\n    features_ts_to_create += [f\"abs_{f}_mean\", f\"abs_{f}_std\"]\n\nfor f in [\"anglez\", \"enmo\", \"X\", \"Y\", \"Z\", \"norm_xyz\", \"norm_xy\", \"norm_xz\", \"norm_yz\", \"light\"]:\n    for per in [\"all\", \"night\", \"evening\"]:\n        features_ts_to_create += [f\"{f}_{per}_{indic_}\" for indic_ in [\"q75\", \"q90\", \"med\", \"mean\", \"std\", \"skew\", \"kurt\"]]\n\nfor f in [\"anglez\", \"enmo\", \"X\", \"Y\", \"Z\", \"norm_xyz\", \"norm_xy\", \"norm_xz\", \"norm_yz\", \"light\"]:\n    for per in [\"all\", \"night\", \"evening\"]:\n        for q in [\"q10\", \"q25\", \"med\", \"q75\", \"q90\", \"q95\"]:\n            features_ts_to_create += [f\"{q}_max_{f}_{per}\", f\"{q}_med_{f}_{per}\"]\n\n#for f in [\"anglez\", \"enmo\", \"X\", \"Y\", \"Z\", \"norm_xyz\", \"norm_xy\", \"norm_xz\", \"norm_yz\", \"light\"]:\nfor f in [\"anglez\", \"enmo\", \"norm_xyz\", \"norm_xy\", \"norm_xz\", \"norm_yz\", \"light\"]:\n    for per in [\"school\", \"afternoon\", \"wevening\", \"weekendday\"]:\n        for q in [\"q10\", \"q25\", \"med\", \"q75\", \"q90\", \"q95\"]:\n            features_ts_to_create += [f\"{q}_max_{f}_{per}\", f\"{q}_med_{f}_{per}\"]\n\nfor f in [\"bedtime\", \"waketime\", \"sleep_hours\"]:\n    for q in [\"mean\", \"std\"]:\n        features_ts_to_create += [f\"{f}_{q}\"]\n    for v in [\"we1\", \"we2\"]:\n        for p in [0, 1]:\n            for q in [\"mean\", \"std\"]:\n                features_ts_to_create += [f\"{q}_{v}_{p}_{f}\"]\n\nfeatures_ts_to_create += ['enmo_X_light_we_mean']\n\nprint(f\"There were first {len(features_ts_to_create)} available time Series features\\n\\n\")\n\nfeatures_ts_to_create = [f for f in features_ts_to_create if f in features_to_create] # Order preserved\n\nprint(f\"But finallly, there are {len(features_ts_to_create)} time Series features to create : \\n{features_ts_to_create} \\n\\n\")\n\ndef func_sleep_wake(x, sleephr, awakeval, wakehr):\n    '''\n    # Fit Function used to fit bedtime and waketime.\n    # Returns a function of x (time in hours) with three regions:\n    #    -------\\________/------------\n    #  awakeval  sleepval    awakeval\n    # sleepval and a ramp transistion time are fixed in the code.\n    # Most useful if the x values go from awake time to awake time,\n    # e.g., 3 pm (x=15 hours) to the next 3 pm (x=15+24 hours).\n    '''\n    sleepval = 0.05\n    ramp = 4.0*(1/6)  # hours\n    out = (sleepval + (awakeval - sleepval) *\n                            (1.0 - np.clip((x - (sleephr - 0*ramp))/ramp, 0.0,1.0) *\n                                 np.clip(((wakehr + ramp/2) - x)/ramp, 0.0,1.0)))\n    return out\n    \ndef process_file(filename, dirname):\n    \n    all_features = []\n    \n    df = pd.read_parquet(os.path.join(dirname, filename, 'part-0.parquet')).drop('step', axis = 1)\n    df[\"hour\"] = (np.trunc(24 * df['time_of_day'] / 86400e9)).astype(np.int8)\n    #df['hour_of_day'] = df['time_of_day'] / 1e9 / 3600   # 0 to 23.9999 \n    #print(24 / 86400e9, 1 /1e9 / 3600 )    # Hour and hour of day are equal\n    \n    # ndays\n    all_features.append((df[\"relative_date_PCIAT\"].max() - df[\"relative_date_PCIAT\"].min()).astype(np.int32))\n\n    # ndays with more than 5000 points\n    p_ = df.groupby(\"relative_date_PCIAT\")[\"enmo\"].count()\n    all_features.append(p_[p_>5000].shape[0])\n    \n    # ndays with max light higher than 1000\n    h_ = df.groupby(\"relative_date_PCIAT\")[\"light\"].max()\n    all_features.append(h_[(h_>1000) & (p_>5000)].shape[0])\n\n    # median max light\n    all_features.append(h_[p_>5000].median())\n    \n    # ndays with max light higher than 2000\n    all_features.append(h_[(h_>2000)& (p_>5000)].shape[0])\n    \n    # ndays with max light lower than 1000\n    all_features.append(h_[(h_<1000) & (p_>5000)].shape[0])\n    \n    # Max light\n    all_features.append(df[\"light\"].max())\n\n    # Right-handed\n    all_features.append(int(df[\"X\"].mean()<0))\n\n    for f in [\"enmo\", \"X\", \"Y\", \"Z\", \"anglez\"]:\n        df[f\"abs_{f}\"] = df[f].abs()\n        h_ = df.groupby(\"relative_date_PCIAT\")[f\"abs_{f}\"].mean()\n        if f\"abs_{f}_mean\" in features_to_create: all_features.append(h_[p_>5000].mean())\n        if f\"abs_{f}_std\" in features_to_create: all_features.append(h_[p_>5000].std())\n\n    # Many features by loop \n    df[\"norm_xyz\"] = np.sqrt(df[\"X\"]**2 + df[\"Y\"]**2 + df[\"Z\"]**2)\n    df[\"norm_xy\"]  = np.sqrt(df[\"X\"]**2 + df[\"Y\"]**2)\n    df[\"norm_xz\"]  = np.sqrt(df[\"X\"]**2 + df[\"Z\"]**2)\n    df[\"norm_yz\"]  = np.sqrt(df[\"Y\"]**2 + df[\"Z\"]**2)\n    \n    # masks\n    mask_worne = (df['non-wear_flag'] == 0)\n    \n    mask_night = df[\"hour\"].isin([22, 23, 0, 1, 2, 3, 4, 5, 6])\n    mask_night &= (df['non-wear_flag'] == 0)\n    \n    mask_evening = df[\"hour\"].isin([16, 17, 18, 19, 20, 21])\n    mask_evening &= (df['non-wear_flag'] == 0)\n    \n    for f in [\"anglez\", \"enmo\", \"X\", \"Y\", \"Z\", \"norm_xyz\", \"norm_xy\", \"norm_xz\", \"norm_yz\", \"light\"]:\n        for mask, per in zip([mask_worne, mask_night, mask_evening], [\"all\", \"night\", \"evening\"]):\n            if f\"{f}_{per}_q75\"  in features_to_create: all_features.append(df.loc[mask, f].quantile(q =.75))\n            if f\"{f}_{per}_q90\"  in features_to_create: all_features.append(df.loc[mask, f].quantile(q =.9))\n            if f\"{f}_{per}_med\"  in features_to_create: all_features.append(df.loc[mask, f].quantile(q =.5))\n            if f\"{f}_{per}_mean\" in features_to_create: all_features.append(df.loc[mask, f].mean())\n            if f\"{f}_{per}_std\"  in features_to_create: all_features.append(df.loc[mask, f].std())\n            if f\"{f}_{per}_skew\" in features_to_create: all_features.append(df.loc[mask, f].skew())\n            if f\"{f}_{per}_kurt\" in features_to_create: all_features.append(df.loc[mask, f].kurt())\n    \n    for f in [\"anglez\", \"enmo\", \"X\", \"Y\", \"Z\", \"norm_xyz\", \"norm_xy\", \"norm_xz\", \"norm_yz\", \"light\"]:\n        for mask, per in zip([mask_worne, mask_night, mask_evening], [\"all\", \"night\", \"evening\"]):\n            h_ = df.loc[mask].groupby(\"relative_date_PCIAT\")[f].agg([\"max\", \"median\"])\n            for q, libq in zip([.1, .25, .5, .75, .9, .95], [\"q10\", \"q25\", \"med\", \"q75\", \"q90\", \"q95\"]):\n                if f\"{libq}_max_{f}_{per}\" in features_to_create: all_features.append(h_[\"max\"][p_>5000].quantile([q]).values[0])\n                if f\"{libq}_med_{f}_{per}\" in features_to_create: all_features.append(h_[\"median\"][p_>5000].quantile([q]).values[0])\n\n    mask_wd = df[\"weekday\"].isin([1, 2, 3, 4, 5])\n    mask_we = df[\"weekday\"].isin([6, 7])\n    \n    mask_school = df[\"hour\"].isin([7, 8, 9, 10, 11, 12, 13, 14])\n    mask_school &= (df['non-wear_flag'] == 0)\n    mask_school &= mask_wd\n    \n    mask_afternoon = df[\"hour\"].isin([15, 16, 17, 18, 19])\n    mask_afternoon &= (df['non-wear_flag'] == 0)\n    mask_afternoon &= mask_wd\n    \n    mask_evening = df[\"hour\"].isin([20, 21])\n    mask_evening &= (df['non-wear_flag'] == 0)\n    mask_evening &= mask_wd\n\n    mask_weday = df[\"hour\"].isin([7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21])\n    mask_weday &= (df['non-wear_flag'] == 0)\n    mask_weday &= mask_we\n    \n#    for f in [\"anglez\", \"enmo\", \"X\", \"Y\", \"Z\", \"norm_xyz\", \"norm_xy\", \"norm_xz\", \"norm_yz\", \"light\"]:\n    for f in [\"anglez\", \"enmo\", \"norm_xyz\", \"norm_xy\", \"norm_xz\", \"norm_yz\", \"light\"]:\n        for mask, per in zip([mask_school, mask_afternoon, mask_evening, mask_weday], [\"school\", \"afternoon\", \"wevening\", \"weekendday\"]):\n            h_ = df.loc[mask].groupby(\"relative_date_PCIAT\")[f].agg([\"max\", \"median\"])\n            for q, libq in zip([.1, .25, .5, .75, .9, .95], [\"q10\", \"q25\", \"med\", \"q75\", \"q90\", \"q95\"]):\n                if f\"{libq}_max_{f}_{per}\" in features_to_create: all_features.append(h_[\"max\"][p_>5000].quantile([q]).values[0])\n                if f\"{libq}_med_{f}_{per}\" in features_to_create: all_features.append(h_[\"median\"][p_>5000].quantile([q]).values[0])\n\n    check_bedtime = False\n    for f in [\"bedtime\", \"waketime\", \"sleep_hours\"]:\n        if f\"{f}_mean\" in features_ts_to_create \\\n        or f\"{f}_std\" in features_ts_to_create \\\n        or f\"mean_we1_0_{f}\" in features_ts_to_create \\\n        or f\"mean_we1_1_{f}\" in features_ts_to_create \\\n        or f\"mean_we2_0_{f}\" in features_ts_to_create \\\n        or f\"mean_we2_1_{f}\" in features_ts_to_create \\\n        or f\"std_we1_0_{f}\" in features_ts_to_create \\\n        or f\"std_we1_1_{f}\" in features_ts_to_create \\\n        or f\"std_we2_0_{f}\" in features_ts_to_create \\\n        or f\"std_we2_1_{f}\" in features_ts_to_create: check_bedtime = True\n        \n    if check_bedtime:\n        # ---89, 90--- Wrist-bedtime, Wrist-sleep_hours\n        # Fit a ----\\___/---- function to determine average bedtime and waketime.\n        # (Could select just weekday sleep/wakes, e.g., Mon 3pm to Fri 3pm.)\n        # Sleep/wake times can vary, but person is almost always awake at 3pm:\n        # Have time go from 15.0 to 24+15 hours - just change the hours, not rearrange\n        enmo995 = df[\"enmo\"].quantile(0.995)\n        df['enmo'] = df[\"enmo\"].clip(0.0, enmo995)\n        df['enmoMed'] = df[\"enmo\"].rolling(121, center = True, min_periods = 1).median()\n        df['enmoMedLog'] = 2.0 + np.log10(0.01 + df['enmoMed'])\n    \n        wrap_hour = df[\"hour\"].copy()\n        sel_wrap = (df[\"hour\"] < 15.0)\n        wrap_hour[sel_wrap] += 24.0\n        wrap_day = df[\"relative_date_PCIAT\"].copy()\n        wrap_day[sel_wrap] += -1\n        wrap_enmo = df['enmoMedLog'] # values and locations don't change\n\n        bedtimes, waketimes, sleep_hours_l, weekdays = [], [], [], []\n#    warnings.filterwarnings(\"ignore\")\n        for nd, d in enumerate(list(wrap_day.unique())):\n#            if nd < 22:\n            if 1 == 1:\n                lower_bounds = [ 16.0,       0.2,   5.0+24.0]\n                upper_bounds = [  5.0+24.0,  3.0,  12.0+24.0]\n                try:\n                    # curve_fit from scipy.optimize\n                    if (wrap_hour.loc[wrap_day == d].count() - wrap_hour.loc[wrap_day == d].isna().sum()) > 10:\n                        bedtime, wake_level, waketime = curve_fit(\n                            func_sleep_wake, wrap_hour.loc[wrap_day == d], wrap_enmo.loc[wrap_day == d], sigma=0.05, nan_policy = 'omit',#'raise',\n                            x_scale = [5,0.5,5], bounds = (lower_bounds, upper_bounds), ftol=.01, xtol=.01, gtol=.01)[0]\n                        sleep_hours = waketime - bedtime\n                        bedtimes.append(bedtime)\n                        waketimes.append(waketime-24)\n                        sleep_hours_l.append(sleep_hours)\n                        weekdays.append(df.loc[wrap_day == d, \"weekday\"].values[0])\n                except:\n                    pass\n        df_ = pd.DataFrame([bedtimes, waketimes, sleep_hours_l, weekdays]).transpose()\n        df_.columns = [\"bedtime\", \"waketime\", \"sleep_hours\", \"weekday\"]\n        df_[\"we1\"] = df_[\"weekday\"].map({1:1, 2:1, 3:1, 4:1, 5:0, 6:0, 7:1})\n        df_[\"we2\"] = df_[\"weekday\"].map({1:1, 2:1, 3:1, 4:1, 5:1, 6:0, 7:0})\n\n        for f in [\"bedtime\", \"waketime\", \"sleep_hours\"]:\n            if f\"{f}_mean\" in features_ts_to_create: all_features.append(df_[f].mean())\n            if f\"{f}_std\" in features_ts_to_create: all_features.append(df_[f].std())\n            for v in [\"we1\", \"we2\"]:\n                for p in [0, 1]:\n                    if f\"mean_{v}_{p}_{f}\" in features_ts_to_create: all_features.append(df_.loc[df_[v]==p, f].mean())\n                    if f\"std_{v}_{p}_{f}\" in features_ts_to_create: all_features.append(df_.loc[df_[v]==p, f].std())\n    \n    \n    \n    if 'enmo_X_light_we_mean' in features_ts_to_create:\n        df[\"enmo_X_light\"] = df[\"enmo\"] * df[\"light\"]\n        all_features.append(df.loc[mask_weday, \"enmo_X_light\"].mean())\n        \n#    display(df_)\n\n    return all_features, filename.split('=')[1]\n    \n        \ndef load_time_series(dirname) -> pd.DataFrame:\n    \n    ids = os.listdir(dirname)\n    if debug: ids = ids[:10]\n    \n    with ThreadPoolExecutor() as executor:\n        results = list(tqdm(executor.map(lambda fname: process_file(fname, dirname), ids), total = len(ids)))\n        \n    stats, idx = zip(*results)\n    \n    cols = [\"ndays\", \"ndays_with_more_than_5000_points\", \"ndays_with_light_higher_than_1000\", \"median_max_light\", \n            \"ndays_with_light_higher_than_2000\", \"ndays_with_light_lower_than_1000\", \"max_light\", \"right_handed\"]\n\n    df = pd.DataFrame(stats, columns = cols + features_ts_to_create, index = pd.Index(idx, name = \"id\"))\n    \n    df[\"ratio_2000\"] = df[\"ndays_with_light_higher_than_2000\"] / df[\"ndays_with_more_than_5000_points\"]\n                \n    return df\n\ntrain_ts = load_time_series(f\"{path}/child-mind-institute-problematic-internet-use/series_train.parquet\")\ntest_ts = load_time_series(f\"{path}/child-mind-institute-problematic-internet-use/series_test.parquet\")\n\n#train_ts.to_csv(\"train_ts.csv\")\n#test_ts.to_csv(\"test_ts.csv\")\n\ndisplay(train_ts.head())\ndisplay(test_ts.head())","metadata":{"papermill":{"duration":105.630433,"end_time":"2024-09-25T20:44:08.050521","exception":false,"start_time":"2024-09-25T20:42:22.420088","status":"completed"},"tags":[],"trusted":true,"_kg_hide-input":true,"execution":{"iopub.status.busy":"2024-12-17T21:01:26.542364Z","iopub.execute_input":"2024-12-17T21:01:26.542739Z","iopub.status.idle":"2024-12-17T21:09:03.059364Z","shell.execute_reply.started":"2024-12-17T21:01:26.542708Z","shell.execute_reply":"2024-12-17T21:09:03.058415Z"},"jupyter":{"source_hidden":true}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## Main dataset & FE\n\nCredit : [M abdullah](https://www.kaggle.com/abdullah0a) for features ideas in https://www.kaggle.com/code/abdullah0a/ensamble-models","metadata":{}},{"cell_type":"code","source":"train = pd.read_csv(f'{path}/child-mind-institute-problematic-internet-use/train.csv', index_col = \"id\")\ntest = pd.read_csv(f'{path}/child-mind-institute-problematic-internet-use/test.csv', index_col = \"id\")\n\nfeatures_season = list(test.dtypes[test.dtypes == object].index)\ntarget = \"sii\"\n\nfor df in [train, test]:\n    df['BMI_Age'] = df['Physical-BMI'] * df['Basic_Demos-Age']\n    df['Internet_Hours_Age'] = df['PreInt_EduHx-computerinternet_hoursday'] * df['Basic_Demos-Age']\n    df['BMI_Internet_Hours'] = df['Physical-BMI'] * df['PreInt_EduHx-computerinternet_hoursday']\n    df['BFP_BMI'] = df['BIA-BIA_Fat'] / df['BIA-BIA_BMI']\n    df['FFMI_BFP'] = df['BIA-BIA_FFMI'] / df['BIA-BIA_Fat']\n    df['FMI_BFP'] = df['BIA-BIA_FMI'] / df['BIA-BIA_Fat']\n    df['LST_TBW'] = df['BIA-BIA_LST'] / df['BIA-BIA_TBW']\n    df['BFP_BMR'] = df['BIA-BIA_Fat'] * df['BIA-BIA_BMR']\n    df['BFP_DEE'] = df['BIA-BIA_Fat'] * df['BIA-BIA_DEE']\n    df['BMR_Weight'] = df['BIA-BIA_BMR'] / df['Physical-Weight']\n    df['DEE_Weight'] = df['BIA-BIA_DEE'] / df['Physical-Weight']\n    df['SMM_Height'] = df['BIA-BIA_SMM'] / df['Physical-Height']\n    df['Muscle_to_Fat'] = df['BIA-BIA_SMM'] / df['BIA-BIA_FMI']\n    df['Hydration_Status'] = df['BIA-BIA_TBW'] / df['Physical-Weight']\n    df['ICW_TBW'] = df['BIA-BIA_ICW'] / df['BIA-BIA_TBW']\n\n# Debugging mode, my public score is so low, that I need to ensure test predictions will be predict as oof predictions\n# search \"good line\" to find other changes\nif 2==1:\n    test = train.dropna(subset=target).copy()\n    test_ts = train_ts.copy()\n    print(test.shape)","metadata":{"papermill":{"duration":0.318676,"end_time":"2024-09-25T20:44:08.42056","exception":false,"start_time":"2024-09-25T20:44:08.101884","status":"completed"},"tags":[],"trusted":true,"_kg_hide-input":true,"execution":{"iopub.status.busy":"2024-12-17T21:28:11.744828Z","iopub.execute_input":"2024-12-17T21:28:11.745227Z","iopub.status.idle":"2024-12-17T21:28:11.822996Z","shell.execute_reply.started":"2024-12-17T21:28:11.745192Z","shell.execute_reply":"2024-12-17T21:28:11.821819Z"},"jupyter":{"source_hidden":true}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Functions\n## Feature importance","metadata":{"papermill":{"duration":0.047674,"end_time":"2024-09-25T20:44:08.545066","exception":false,"start_time":"2024-09-25T20:44:08.497392","status":"completed"},"tags":[]}},{"cell_type":"code","source":"class FeatureImportance():\n    \n    def __init__(self, n_permutations = None, lib_metric = None, digit = 4):\n        \n        self.n_permutations = n_permutations\n        self.lib_metric = lib_metric\n        self.differences = {}\n        self.digit = digit\n        \n        self.metric_max = [\"Accuracy\", \"AUC\", \"Recall@20\", \"cohen_kappa_score\", \"auc\", \"f1\", \"R2\", \"QWK\"]\n        self.metric_min = [\"Logloss\", \"rmse\", \"RMSE\", \"mae\", \"Balanced log loss\", \"RMSLE\"]\n        if lib_metric not in self.metric_max + self.metric_min:\n            raise ValueError(f\"Unknown lib_metric : {lib_metric}.\")\n        \n    def append(self, feature, metric_value, metric_ref):\n        \n        if feature not in self.differences.keys():\n            self.differences[feature] = []\n            \n        self.differences[feature].append(metric_value - metric_ref)\n        \n    def build_dataframe(self):\n        \n        if self.differences == {}: return\n        \n        df_diff = pd.DataFrame(self.differences).transpose()\n        self.total_n_permutations = df_diff.shape[1]\n        \n        df_diff[\"mean\"] = df_diff.mean(axis=1)\n        df_diff.sort_values(\"mean\", ascending = self.lib_metric in self.metric_max, inplace = True)\n        \n        for i in range(self.total_n_permutations):\n            if self.lib_metric in self.metric_max: \n                useless_ind = pd.Series((df_diff[i] > 0).values, index = df_diff.index, name = f\"useless_{i}\")\n            else:\n                useless_ind = pd.Series((df_diff[i] < 0).values, index = df_diff.index, name = f\"useless_{i}\")\n            df_diff = pd.concat([df_diff, useless_ind], axis = 1)\n        df_diff[\"useless\"] = df_diff[[f\"useless_{i}\" for i in range(self.total_n_permutations)]].sum(axis=1)\n        \n        return df_diff\n        \n    def plot(self, title=\"\", remove_from_output = []):\n        \n        if self.differences == {}: return\n        \n        df_diff = self.build_dataframe()\n        if len(remove_from_output)>0:\n            df_diff = df_diff.loc[~df_diff.index.isin(remove_from_output)]\n        \n        fk = next(iter(self.differences))\n        len_fk = len(self.differences[fk])\n        \n        fig, ax = plt.subplots(1, 2, figsize = (15, int(df_diff.shape[0] * 1/2)))\n        labels = [f\"{l} ({v:,.0f})\" for l, v in zip(df_diff.index[::-1], df_diff[\"useless\"][::-1])]\n        ax[0].boxplot(df_diff[[c for c in range(len_fk)]][::-1].T,\n                          vert=False, labels=labels)\n        ax[0].axvline(x = 0, color = 'green')\n        ax[0].set_xlabel(f\"Dist. {self.lib_metric} difference after permutation of values (n times useless / {self.total_n_permutations})\")\n            \n        sns.barplot(x = \"mean\", y = df_diff.index, data = df_diff, ax = ax[1], color = \"green\")\n        ax[1].bar_label(ax[1].containers[0], fmt = f'%.{self.digit}f', padding = 2)\n        ax[1].set_xlabel(f\"Mean {self.lib_metric} difference after permutation of values\")\n        ax[1].set_yticklabels([])\n        ax[1].set_yticks([])\n        \n        for ax_ in ax:\n            ax_.set_title(title)\n            ax_.xaxis.set_ticks_position(\"top\")\n            ax_.xaxis.set_label_position('top')\n            ax_.spines[[\"right\", \"bottom\"]].set_visible(False)\n        ax[1].spines[[\"left\"]].set_visible(False)","metadata":{"_kg_hide-input":true,"papermill":{"duration":0.07501,"end_time":"2024-09-25T20:44:08.668881","exception":false,"start_time":"2024-09-25T20:44:08.593871","status":"completed"},"tags":[],"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T21:28:12.075478Z","iopub.execute_input":"2024-12-17T21:28:12.076227Z","iopub.status.idle":"2024-12-17T21:28:12.093136Z","shell.execute_reply.started":"2024-12-17T21:28:12.076188Z","shell.execute_reply":"2024-12-17T21:28:12.091991Z"},"jupyter":{"source_hidden":true}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## More parameters to GBM","metadata":{}},{"cell_type":"code","source":"def complete_params(regressor, params):\n    \n    model = regressor(**params)\n\n    new_params = deepcopy(params)\n    \n    if isinstance(model, LGBMRegressor):\n        new_params.update({\n            \"device\"                : \"gpu\" if cuda else \"cpu\",\n            \"metric\"                : \"rmse\",\n            \"random_state\"          : seed,\n            \"n_estimators\"          : 2000, 'early_stopping_round' : 100,\n            \"n_jobs\"                : -1,\n            \"verbose\"               : -1,\n        })\n    if isinstance(model, CatBoostRegressor):\n        new_params.update({\n            'task_type'           : \"GPU\" if cuda else \"CPU\",\n            'n_estimators'        : 2000, 'od_wait': 100, \n            'loss_function'       : \"RMSE\", #'MultiClass', \n            'eval_metric'         : \"RMSE\",#\"Accuracy\",  # à partir de cb3\n#            'score_function'      : \"NewtonL2\",\n            'random_state'        : seed,\n            'use_best_model'      : True,\n            'grow_policy'         : \"SymmetricTree\", # default for catboost\n        })\n    if isinstance(model, XGBRegressor):\n        new_params.update({\n            'booster'               : 'gbtree',\n            \"tree_method\"           : \"hist\",\n            'objective'             : 'reg:squarederror', #'binary:logistic'\n            'eval_metric'           : \"rmse\", #\"auc\",\n            \"device\"                : \"cuda\" if cuda else \"cpu\",\n            'verbosity'             : 0,\n            \"random_state\"          : seed,\n            'n_estimators'          : 3000, 'early_stopping_rounds' : 100,\n            'enable_categorical'    : True,\n        })\n    return new_params","metadata":{"trusted":true,"_kg_hide-input":true,"jupyter":{"source_hidden":true},"execution":{"iopub.status.busy":"2024-12-17T21:28:12.391345Z","iopub.execute_input":"2024-12-17T21:28:12.391763Z","iopub.status.idle":"2024-12-17T21:28:12.399952Z","shell.execute_reply.started":"2024-12-17T21:28:12.391725Z","shell.execute_reply":"2024-12-17T21:28:12.398874Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## Stratification","metadata":{}},{"cell_type":"code","source":"split_class = \"class\"\n\ndef create_class_to_stratify(verbose = False, split_class = split_class):\n    '''\n    Stratification on target and the three most important features\n    ''' \n    \n    # Most important features or columns\n    df = train[[target, \"SDS-SDS_Total_Raw\", \"Basic_Demos-Age\", \"PreInt_EduHx-computerinternet_hoursday\"]].dropna(subset=target).copy()\n    df[\"ndays\"] = train_ts[\"ndays\"]\n\n    # 33th, 66th percentiles of numeric columns\n    for f in [\"SDS-SDS_Total_Raw\", \"Basic_Demos-Age\", \"ndays\"]:\n        df[f\"q_{f}\"] = pd.qcut(df[f], 3).astype(str)\n        if verbose: \n            display(df[f\"q_{f}\"].value_counts(dropna=False))\n    if verbose:\n        display(pd.crosstab(df[target], df[\"q_ndays\"], dropna = False))\n\n    # Create class for stratification\n    class_enc = OrdinalEncoder()\n    df = class_enc.fit_transform(df[[\n        target, \"q_SDS-SDS_Total_Raw\", \"q_Basic_Demos-Age\", \"PreInt_EduHx-computerinternet_hoursday\", \"q_ndays\"\n    ]].astype(str).fillna(-1))\n\n    df[split_class] = 0\n    for i, f in enumerate([target, \"q_SDS-SDS_Total_Raw\", \"q_Basic_Demos-Age\", \"PreInt_EduHx-computerinternet_hoursday\", \"q_ndays\"]):\n        df[split_class] += df[f]**(i*5)\n\n    # Check \n    mask1 = df[split_class].value_counts(dropna=False)\n    mask = df[split_class].isin(mask1.loc[mask1 < n_splits].index)\n    df.loc[mask, split_class] = -1\n    if verbose: display(df[split_class].value_counts(dropna=False))\n\n    return df[split_class]\n\nif split_class in train.columns:\n    train.drop([split_class], axis=1, inplace= True)\n\ntrain = pd.concat([train, create_class_to_stratify()], axis=1)\n#train[\"class\"].value_counts(dropna=False)","metadata":{"_kg_hide-output":true,"trusted":true,"_kg_hide-input":true,"jupyter":{"source_hidden":true},"execution":{"iopub.status.busy":"2024-12-17T21:28:12.734672Z","iopub.execute_input":"2024-12-17T21:28:12.735065Z","iopub.status.idle":"2024-12-17T21:28:12.787107Z","shell.execute_reply.started":"2024-12-17T21:28:12.735028Z","shell.execute_reply":"2024-12-17T21:28:12.786008Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## Threshold\nCredit : [AmbrosM](https://www.kaggle.com/ambrosm) for threshold optimizer in https://www.kaggle.com/code/ambrosm/piu-eda-which-makes-sense","metadata":{"papermill":{"duration":0.048555,"end_time":"2024-09-25T20:44:08.768064","exception":false,"start_time":"2024-09-25T20:44:08.719509","status":"completed"},"tags":[]}},{"cell_type":"code","source":"def round_with_thresholds(raw_preds, thresholds):\n    \"\"\"Round the raw predictions using specified thresholds\n    \n    Parameters\n    ----------\n    raw_preds: raw predictions of the regressor, array of n_samples float values\n    thresholds: 3-element float array\n    \n    Returns\n    -------\n    rounded_preds: rounded predictions, array of n_samples int values in range 0..3\n    \"\"\"\n    return np.where(raw_preds < thresholds[0], 0,\n                    np.where(raw_preds < thresholds[1], 1,\n                             np.where(raw_preds < thresholds[2], 2, 3)))\n\n\ndef fun(thresholds, y_true, raw_preds):\n    \"\"\"Function to be minimized: negative quadratic kappa score\n    \n    Parameters:\n    thresholds: ndarray of shape (3, )\n    y_true: ndarray of shape (n_samples, )\n    raw_preds: ndarray of shape (n_samples, )\n    \n    Returns:\n    negative quadratic kappa score for the predictions rounded at the specified thresholds\n    \"\"\"\n    rounded_preds = round_with_thresholds(raw_preds, thresholds)\n    return - cohen_kappa_score(y_true, rounded_preds, weights='quadratic')","metadata":{"trusted":true,"jupyter":{"source_hidden":true},"_kg_hide-input":true,"execution":{"iopub.status.busy":"2024-12-17T21:28:13.070077Z","iopub.execute_input":"2024-12-17T21:28:13.070466Z","iopub.status.idle":"2024-12-17T21:28:13.077536Z","shell.execute_reply.started":"2024-12-17T21:28:13.070432Z","shell.execute_reply":"2024-12-17T21:28:13.076386Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## To categorical encoder","metadata":{}},{"cell_type":"code","source":"class to_categorical(BaseEstimator, TransformerMixin):\n    \n    def __init__(self, min_frequency = 10, unknown_value = -666, nan_value = -667):\n        \n        self.verbose = True\n        self.min_frequency = min_frequency\n        self.unknown_value = unknown_value\n        self.nan_value = nan_value\n            \n    def fit(self, X, y = None):\n        \n        df = X.copy()\n        self.all_categ = {}\n        self.transf = {}\n        \n        for f in df.columns:\n            \n            self.transf[f] = []\n            vc = df[f].value_counts()\n            for v in vc[vc < self.min_frequency].index:\n                df.loc[df[f] == v, f] = self.unknown_value\n                self.transf[f].append(v)\n            \n            self.all_categ[f] = df[f].astype(\"category\").cat\n            \n            if self.unknown_value not in list(self.all_categ[f].categories):\n                self.all_categ[f] = self.all_categ[f].add_categories([self.unknown_value]).cat\n            self.all_categ[f] = self.all_categ[f].add_categories([self.nan_value]).cat\n\n        return self\n    \n    \n    def transform(self, X, y = None):\n        \n        df = X.copy()\n        \n        for f in df.columns:\n            \n            for v in self.transf[f]:\n                df.loc[df[f] == v, f] = self.unknown_value\n            \n            df[f] = pd.Categorical(df[f], categories = self.all_categ[f].categories)\n\n            df.loc[df[f].isna(), f] = self.nan_value\n            \n        return df","metadata":{"_kg_hide-input":true,"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T21:28:13.447706Z","iopub.execute_input":"2024-12-17T21:28:13.448589Z","iopub.status.idle":"2024-12-17T21:28:13.459415Z","shell.execute_reply.started":"2024-12-17T21:28:13.448543Z","shell.execute_reply":"2024-12-17T21:28:13.458188Z"},"jupyter":{"source_hidden":true}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"class to_numpy(BaseEstimator, TransformerMixin):\n    \n    def __init__(self):\n        self.verbose = True\n            \n    def fit(self, X, y = None):\n        return self\n    \n    def transform(self, X, y = None):\n        return X.values","metadata":{"_kg_hide-input":true,"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T21:28:13.750614Z","iopub.execute_input":"2024-12-17T21:28:13.751020Z","iopub.status.idle":"2024-12-17T21:28:13.756925Z","shell.execute_reply.started":"2024-12-17T21:28:13.750988Z","shell.execute_reply":"2024-12-17T21:28:13.755867Z"},"jupyter":{"source_hidden":true}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## CV Loop function","metadata":{}},{"cell_type":"code","source":"def do_cv(train = train.dropna(subset=target), target = target, test = test, other_data = None, features = [], digit = 3,\n          enc = ColumnTransformer(transformers = [], remainder = \"passthrough\", verbose_feature_names_out = False), \n          my_model = None, params = {}, do_feat_imp = False, metric_label = \"QWK\", log_eval_for_gbms=2000,\n          folds = RepeatedStratifiedKFold(n_repeats = n_repeats, n_splits = n_splits, random_state = seed), \n          split_class = target, target_mapper = None, outlier_pipeline = None, ndays_min = None, label = None, \n          x0 = [0.5, 1.5, 2.5],\n         ):\n\n    random_state, n_repeats = folds.__dict__[\"random_state\"], folds.__dict__[\"n_repeats\"]\n    \n    feat_imp = FeatureImportance(n_permutations = 2, lib_metric = metric_label)\n    oof_scores, val_scores, trn_scores, best_iteration, history = [], [], [], [], {}\n    oofs, preds = pd.Series(0.0, name = target, index = train.index), pd.Series(0.0, name = target, index = test.index)\n    df_max_bins = pd.DataFrame()\n    \n    for fold, (trn_idx, val_idx) in enumerate(folds.split(train, train[split_class])):\n\n        # Training samples\n        if other_data is not None:      # if some augmentations...\n            X_trn = pd.concat([train.iloc[trn_idx], other_data], axis=0)\n            y_trn = pd.concat([train[target].iloc[trn_idx], other_data[target]], axis=0)\n        else:\n            X_trn, y_trn = train.iloc[trn_idx], train[target].iloc[trn_idx]\n            \n        # Validation samples\n        X_val, y_val = train.iloc[val_idx], train[target].iloc[val_idx]\n\n        if ndays_min is not None:\n            X_trn = X_trn.merge(train_ts.loc[train_ts[\"ndays_with_more_than_5000_points\"] >= ndays_min], \n                                how = \"left\", left_index = True, right_index = True)[features]\n            X_val = X_val.merge(train_ts.loc[train_ts[\"ndays_with_more_than_5000_points\"] >= ndays_min], \n                                how = \"left\", left_index = True, right_index = True)[features]\n            X_test = test.merge(test_ts.loc[test_ts[\"ndays_with_more_than_5000_points\"] >= ndays_min], \n                                how = \"left\", left_index = True, right_index = True)[features]\n        else:\n            X_trn = X_trn.merge(train_ts, how = \"left\", left_index = True, right_index = True)[features]\n            X_val = X_val.merge(train_ts, how = \"left\", left_index = True, right_index = True)[features]\n            X_test = test.merge(test_ts, how = \"left\", left_index = True, right_index = True)[features]\n            \n        if outlier_pipeline is not None: # to drop outliers\n            outliers = outlier_pipeline.fit_predict(X_trn)\n            X_trn = X_trn[outliers == 1] ; y_trn = y_trn[outliers == 1]\n            \n        y_trn_, y_val_ = y_trn, y_val   \n        if target_mapper is not None:\n            y_trn_, y_val_ = y_trn.map(target_mapper), y_val.map(target_mapper)\n            \n        if (fold) % n_splits == 0:\n            oof_pred, pred = pd.Series(0.0, name = target, index = train.index), pd.Series(0.0, name = target, index = test.index)\n        \n        # Training\n        model = my_model(**params)\n        if \"early_stopping_round\" in params.keys() or (\"early_stopping_rounds\" in params.keys()) or ('od_wait' in params.keys()):\n            if isinstance(model, CatBoostRegressor):\n                model.fit(enc.fit_transform(X_trn), y_trn, eval_set=[(enc.transform(X_val), y_val)], verbose = log_eval_for_gbms)\n                best_iteration.append(model.get_best_iteration())\n            elif isinstance(model, XGBRegressor):\n                model.fit(enc.fit_transform(X_trn), y_trn, eval_set=[(enc.transform(X_val), y_val)], verbose = log_eval_for_gbms)\n                best_iteration.append(model.best_iteration)\n                history = model.evals_result()\n            elif isinstance(model, LGBMRegressor):\n                model.fit(enc.fit_transform(X_trn), y_trn_, eval_set=[(enc.transform(X_val), y_val_)], \n                          callbacks = [early_stopping(params[\"early_stopping_round\"], verbose = False), log_evaluation(log_eval_for_gbms)])\n                best_iteration.append(model.best_iteration_)\n            else: print(\"Error\")\n        else:\n            if isinstance(model, CatBoostRegressor) or isinstance(model, XGBRegressor) or isinstance(model, LGBMRegressor):\n                model.fit(enc.fit_transform(X_trn), y_trn, eval_set=[(enc.transform(X_trn), y_trn), (enc.transform(X_val), y_val)], verbose = log_eval_for_gbms)\n                best_iteration.append(params[\"n_estimators\"])\n            else:\n                model.fit(enc.fit_transform(X_trn), y_trn)\n            if isinstance(model, XGBRegressor):\n                history[f\"fold{fold+1}\"] = model.evals_result()\n\n        # Score & prediction\n        y_trn_pred = model.predict(enc.transform(X_trn))\n        y_val_pred = model.predict(enc.transform(X_val))\n        y_pred     = model.predict(enc.transform(X_test))\n        pred += y_pred / n_splits              # Good line \n#        pred.iloc[val_idx] += y_pred[val_idx]  # Bad line : just to test predictions with train dataset\n        oof_pred.iloc[val_idx] += y_val_pred\n        oofs.iloc[val_idx] += oof_pred.iloc[val_idx] / n_repeats\n        \n        min_bins = minimize(fun, x0 = x0, args = (y_trn, y_trn_pred), method = 'Nelder-Mead')\n        df_max_bins = pd.concat([df_max_bins, pd.DataFrame(min_bins.x).transpose()], axis=0)\n        trn_scores.append(score_(y_trn, round_with_thresholds(y_trn_pred, min_bins.x)))\n        val_scores.append(score_(y_val, round_with_thresholds(y_val_pred, min_bins.x)))\n        \n        print(f\"    Fold {fold + 1:2}\", end='')\n        if outlier_pipeline is not None:\n            print(f\" (trained without {outliers[outliers==-1].shape[0]} outliers)\", end='')\n        print(f\" : {metric_label} {val_scores[-1:][0]:.{digit}f}\", end='')\n        print(f\" | in train {trn_scores[-1:][0]:.{digit}f}\", end='')\n        print(f' | Overfitting {trn_scores[-1:][0] - val_scores[-1:][0]:.{digit}f}', end='')\n        if (len(best_iteration) > 0) and (best_iteration[-1:][0] > 0):\n            print(f' | Best iteration {best_iteration[-1:][0]}')\n        else:\n            print('')\n            \n        # Mean of scores after n_splits\n        if (fold + 1) % n_splits == 0:\n            \n            oof_label = pd.Series(\n                round_with_thresholds(oof_pred, df_max_bins[-n_splits::].mean(axis=0).values), index = train.index, name = target)\n            oof_scores.append(score_(train[target], oof_label))\n            oof_label_ = pd.Series(\n                round_with_thresholds(oofs * n_repeats / ((fold + 1 ) // n_splits), df_max_bins.mean(axis=0).values), \n                index = train.index, name = target\n            ) \n            pred_label = pd.Series(\n                round_with_thresholds(pred, df_max_bins[-n_splits::].mean(axis=0).values), index = test.index, name = target)\n            \n            print(f'  {bold}Repeat {(fold + 1 ) // n_splits:2} : OOF {metric_label} {bold_blue}{oof_scores[-1]:.{digit}f}{end} | ', end='')\n            print(f'Mean {metric_label} {np.mean(val_scores[-n_splits:]):.{digit}f}({bold+Red}±{np.std(val_scores[-n_splits:]):.{digit}f}){end}) | ', end='')\n            print(f'Cum. {score_(train[target], oof_label_):.{digit}f}{end}')\n            \n            if label is not None:\n#                df_max_bins[-n_splits::].mean(axis=0).to_csv(f\"{output_path}/thresholds_{label}_r{(fold + 1) // n_splits}.csv\")\n#                oof_pred.to_csv(f\"{output_path}/oof_{label}_r{(fold + 1) // n_splits}.csv\")\n#                pred.to_csv(f\"{output_path}/pred_{label}_r{(fold + 1) // n_splits}.csv\")\n                oof_label.to_csv(f\"{output_path}/oof_label_{label}_r{(fold + 1) // n_splits}.csv\")\n                pred_label.to_csv(f\"{output_path}/pred_label_{label}_r{(fold + 1) // n_splits}.csv\")\n            \n            preds += pred / n_repeats\n        \n        # Feature importance\n        if do_feat_imp:\n            for feat in X_val.columns:\n                for i in range(feat_imp.n_permutations):\n                    df = X_val.copy()\n                    # Permutation only on samples with no NAN on this feature :\n                    mask = ~df[feat].isna()\n                    df.loc[mask, feat] = pd.Series(np.random.permutation(df.loc[mask, feat]), index=df.loc[mask].index).astype(X_val[feat].dtypes)\n                    y_p = model.predict(enc.transform(df))\n                    sc = score_(y_val, round_with_thresholds(y_p, min_bins.x))\n                    feat_imp.append(feat, sc, val_scores[-1:][0])\n    # end of loop\n\n    oofs_label = pd.Series(round_with_thresholds(oofs, df_max_bins.mean(axis=0).values), index = train.index, name = target)\n    preds_label = pd.Series(round_with_thresholds(preds, df_max_bins.mean(axis=0).values), index = test.index, name = target)\n    oof_score = score_(train[target], oofs_label)\n\n    if label is not None:\n        \n#        oofs.to_csv(f\"{output_path}/oof_{label}.csv\")\n#        preds.to_csv(f\"{output_path}/pred_{label}.csv\")\n#        df_max_bins.mean(axis=0).to_csv(f\"{output_path}/thresholds_{label}.csv\")\n#        oofs_label.to_csv(f\"{output_path}/oof_label_{label}.csv\")\n#        preds_label.to_csv(f\"{output_path}/pred_label_{label}.csv\")\n\n        feat_imp.plot(title = f\"{label} - {oof_score:.{digit}f}\")\n    else:\n        feat_imp.plot(title = f\"{oof_score:.{digit}f}\")\n        \n    print()\n    if label is not None:\n        print(f'{bold}Cumulative OOF {metric_label} for {label:10} : {bold_blue}{oof_score:.{digit}f}{end} | ', end='')\n    else:\n        print(f'{bold}Cumulative OOF {metric_label} : {bold_blue}{oof_score:.{digit}f}{end} | ', end='')\n    print(f'Mean over OOF : {bold_blue}{np.mean(oof_scores):.{digit}f}{end}', end='')\n    print(f'({bold+Red}±{np.std(oof_scores):.{digit}f}){end} | ', end='')\n    print(f'Mean over folds : {bold_blue}{np.mean(val_scores):.{digit}f}{end}({bold+Red}±{np.std(val_scores):.{digit}f}){end}) | ', end='')\n    print(f'Overfitting : {np.mean(trn_scores)-np.mean(val_scores):.{digit}f}')\n    if len(best_iteration)>0:\n        print(f'Mean Best iteration {int(np.mean(best_iteration))}(±{int(np.std(best_iteration))}) and median {np.median(best_iteration)}')\n        \n    res = {'history': history, \"mean_score\" : val_scores, \"feat_imp\": feat_imp, \"oof_scores\":oof_scores, \"oof_score\":oof_score, \n           \"trn_score\":trn_scores}\n    if len(best_iteration) > 0:\n        res.update({\"mean_it\": int(np.mean(best_iteration)), \"median_it\": np.median(best_iteration), \"best_iteration\":best_iteration})\n    return res","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T21:28:14.741234Z","iopub.execute_input":"2024-12-17T21:28:14.741637Z","iopub.status.idle":"2024-12-17T21:28:14.781001Z","shell.execute_reply.started":"2024-12-17T21:28:14.741601Z","shell.execute_reply":"2024-12-17T21:28:14.780143Z"},"_kg_hide-input":true,"jupyter":{"source_hidden":true}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Training","metadata":{}},{"cell_type":"code","source":"all_results = {}\n\nfor i, (label, p) in enumerate(all_params.items()):\n\n    cond = i==0\n    cond = label in [\"xgb26\"]\n    cond = 1==1\n\n    if cond:\n        \n        print('\\n\\n\\n', label)\n\n        new_model = {\n            'label'           : label,\n            'my_model'        : p[\"regressor\"],\n            'params'          : complete_params(p[\"regressor\"], p[\"params\"]),\n            'features'        : p[\"features\"],\n            'enc'             : ColumnTransformer(transformers = [\n                (\"categorical\", to_categorical(), [f for f in p[\"features\"] if f in features_season])\n            ], remainder = \"passthrough\", verbose_feature_names_out = False),\n            \"ndays_min\"       : p[\"ndays_min\"],\n            'target_mapper'   : p[\"target_mapper\"],\n            'folds'           : RepeatedStratifiedKFold(n_repeats = n_repeats, n_splits = n_splits, random_state = seed), \n            'split_class'     : \"class\",\n            'do_feat_imp'     : do_feat_imp,\n        }\n\n        # A pipeline to detect outlier and drop them from training, just for training, not for validation\n        #   Why : \n        #      - dataset is noisy, remove outlier can reduce noise\n        #      - before using IsolationForest, we need to fill nan values for IsolationForest\n        #      - Nan values are replaced, but only for outliers detection\n        #      - Nan values will not be replaced during traning of GBM\n        if \"isolation\" in list(p.keys()):\n            new_model[\"outlier_pipeline\"] = make_pipeline( \n                ColumnTransformer(transformers = [\n                    (\"imputer_categorical\", make_pipeline(SimpleImputer(strategy = \"most_frequent\"), OrdinalEncoder()), [f for f in p[\"features\"] if f in features_season]),\n                    (\"imputer_numerical\",   SimpleImputer(strategy = \"median\"), [f for f in p[\"features\"] if f not in features_season])\n                ], remainder = \"drop\", verbose_feature_names_out = False), # Keep only features\n                to_numpy(),\n                IsolationForest(n_estimators = 200, contamination = p[\"isolation\"], n_jobs = -1, random_state = seed)\n            )\n            \n        all_results[label] = do_cv(**new_model)","metadata":{"trusted":true,"_kg_hide-output":true,"execution":{"iopub.status.busy":"2024-12-17T21:28:15.985468Z","iopub.execute_input":"2024-12-17T21:28:15.985924Z","iopub.status.idle":"2024-12-17T21:32:22.796999Z","shell.execute_reply.started":"2024-12-17T21:28:15.985889Z","shell.execute_reply":"2024-12-17T21:32:22.795946Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Main results of trainings","metadata":{}},{"cell_type":"code","source":"metric_label, digit = \"QWK\", 3\n\nfor label, res in all_results.items():\n\n    oof_score = res[\"oof_score\"]\n    oof_scores = res[\"oof_scores\"]\n    val_scores = res[\"mean_score\"]\n    trn_scores = res[\"trn_score\"]\n    feat_imp = res[\"feat_imp\"]\n    if \"best_iteration\" in list(res.keys()):\n        best_iteration = res[\"best_iteration\"]\n\n    feat_imp.plot(title = f\"{label} - {oof_score:.{digit}f}\")\n        \n    print(f'{bold}Cumulative OOF {metric_label} for {label:10} : {bold_blue}{oof_score:.{digit}f}{end} | ', end='')\n    print(f'Mean over OOF : {bold_blue}{np.mean(oof_scores):.{digit}f}{end}', end='')\n    print(f'({bold+Red}±{np.std(oof_scores):.{digit}f}){end} | ', end='')\n    print(f'Mean over folds : {bold_blue}{np.mean(val_scores):.{digit}f}{end}({bold+Red}±{np.std(val_scores):.{digit}f}){end}) | ', end='')\n    print(f'Overfitting : {np.mean(trn_scores)-np.mean(val_scores):.{digit}f}')\n    if \"mean_it\" in list(res.keys()):\n        print(f'Mean Best iteration {np.mean(best_iteration):.0f}(±{int(np.std(best_iteration))}) and median {np.median(best_iteration):.0f}')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T21:22:40.762370Z","iopub.execute_input":"2024-12-17T21:22:40.763406Z","iopub.status.idle":"2024-12-17T21:22:44.306300Z","shell.execute_reply.started":"2024-12-17T21:22:40.763364Z","shell.execute_reply":"2024-12-17T21:22:44.305218Z"},"_kg_hide-input":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Partial dependency plots\n## Functions","metadata":{}},{"cell_type":"code","source":"def desc_categorical_features(\n    features, train = train, sharey = False, coutinous_threshold = 20,\n    suptitle_fontsize = \"x-large\", suptitle_weight = \"bold\", \n    title_fontsize = \"large\", title_weight = \"normal\"):\n    \n    for f in features:\n\n        fig = plt.figure(figsize = (15, train[f].nunique()), constrained_layout = True)\n        fig.suptitle(f, fontweight = suptitle_weight, fontsize = suptitle_fontsize)\n        gs = fig.add_gridspec(nrows = 1, ncols = 4)\n        ax = [fig.add_subplot(gs[:, :1]), fig.add_subplot(gs[:, 1:])]\n    \n        train[f].value_counts().sort_index(ascending = False).plot.barh(ax = ax[0], sharey = sharey).set(title = \"N obs in Train\")\n        ax[0].bar_label(ax[0].containers[0], fmt = \"{:,.0f}\", padding = 2)\n    \n        if train[target].nunique() == 2:\n    \n            (train.loc[train[target] == 1, f].value_counts(dropna = False) / train[f].value_counts(dropna = False)).sort_index(ascending = False).plot.barh(ax = ax[1], sharey = sharey)\n            ax[1].set(title = f\"% of {target} == 1 in Train\")\n            ax[1].bar_label(ax[1].containers[0], fmt = \"{:.1%}\", padding = 2)\n        \n        \n        elif train[target].nunique() <= coutinous_threshold:\n            \n            cm = pd.crosstab(train[f], train[target])\n            sns.heatmap(\n                data = (cm.transpose() / cm.sum(axis = 1)).transpose(), ax = ax[1],\n                cmap = sns.dark_palette(\"#69d\", reverse = True, as_cmap = True),\n                cbar = False, lw = 0.25, annot = True, fmt = '.0%'\n            ).set(title = f'{f} per {target} (sum of each row = 100%)')\n    \n        else:\n            \n            sns.boxplot(data = train, x = target, y = f, orient = \"h\", ax=ax[1]).set(title = f)\n            ax[1].set_ylabel(target)\n            \n        for i in range(2):\n            ax[i].title.set(fontsize = title_fontsize, fontweight = title_weight)\n            ax[i].spines[[\"right\", \"bottom\"]].set_visible(False)\n            ax[i].xaxis.set_ticks_position(\"top\")\n            ax[i].set_xlabel(\"\")\n            \ndef desc_continous_features(\n    features, train = train, target = target, target_continous_threshold = 10,\n    suptitle_fontsize = \"x-large\", suptitle_weight = \"bold\", \n    title_fontsize = \"large\", title_weight = \"normal\"):\n\n    warnings.filterwarnings(\"ignore\", category=FutureWarning)\n\n    for f in features:\n        \n        fig = plt.figure(figsize = (15, 5), constrained_layout = True)\n        fig.suptitle(f, fontweight = suptitle_weight, fontsize = suptitle_fontsize)\n        gs = fig.add_gridspec(nrows = 1, ncols = 4)\n        ax = [fig.add_subplot(gs[:, :1]), fig.add_subplot(gs[:, 1:])]\n        \n        sns.histplot(train[f], ax = ax[0]).set(title = f)\n        \n        if train[target].nunique() < target_continous_threshold:\n    \n            sns.boxplot(data = train, x = target, y = f, ax=ax[1]).set(title = target)\n            ax[1].set_xlabel(\"\")\n        \n        else:\n            \n            ax[1].hist2d(x = train[f], y = train[target], cmap = \"Blues\")\n            ax[1].set_xlabel(f)\n            ax[1].set_ylabel(target)\n        \n        for i in range(2):\n            ax[i].title.set(fontsize = title_fontsize, fontweight = title_weight)\n            ax[i].spines[[\"right\", \"bottom\"]].set_visible(False)\n            ax[i].xaxis.set_ticks_position(\"top\")","metadata":{"trusted":true,"_kg_hide-input":true,"execution":{"iopub.status.busy":"2024-12-17T19:31:13.072927Z","iopub.execute_input":"2024-12-17T19:31:13.073327Z","iopub.status.idle":"2024-12-17T19:31:13.091684Z","shell.execute_reply.started":"2024-12-17T19:31:13.073292Z","shell.execute_reply":"2024-12-17T19:31:13.090554Z"},"jupyter":{"source_hidden":true}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## Part 1 : from the main dataset","metadata":{}},{"cell_type":"code","source":"desc_categorical_features([\"Basic_Demos-Sex\"])","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T19:31:13.381566Z","iopub.execute_input":"2024-12-17T19:31:13.381982Z","iopub.status.idle":"2024-12-17T19:31:14.667744Z","shell.execute_reply.started":"2024-12-17T19:31:13.381952Z","shell.execute_reply":"2024-12-17T19:31:14.666343Z"},"_kg_hide-input":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"feats = [ f for f in features_to_create if f in list(test.columns) and f not in [\"Basic_Demos-Sex\"]]\ndesc_continous_features(feats, train = train)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T19:31:17.376665Z","iopub.execute_input":"2024-12-17T19:31:17.377615Z","iopub.status.idle":"2024-12-17T19:31:30.929756Z","shell.execute_reply.started":"2024-12-17T19:31:17.377579Z","shell.execute_reply":"2024-12-17T19:31:30.928601Z"},"_kg_hide-input":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# remember that \n#df['DEE_Weight'] = df['BIA-BIA_DEE'] / df['Physical-Weight']\n#df['LST_TBW'] = df['BIA-BIA_LST'] / df['BIA-BIA_TBW']\n\nt_ = train.copy()\nt_[\"DEE_Weight_clipped\"] = t_[\"DEE_Weight\"].clip(0, 50)\nt_[\"LST_TBW_clipped\"] = t_[\"LST_TBW\"].clip(1.1, 1.5)\ndesc_continous_features([\"DEE_Weight_clipped\", \"LST_TBW_clipped\"], train = t_)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T19:31:30.931899Z","iopub.execute_input":"2024-12-17T19:31:30.932367Z","iopub.status.idle":"2024-12-17T19:31:34.182310Z","shell.execute_reply.started":"2024-12-17T19:31:30.932320Z","shell.execute_reply":"2024-12-17T19:31:34.181063Z"},"_kg_hide-input":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## Part 2 : features from actigraphy data","metadata":{}},{"cell_type":"code","source":"t_ = train.merge(train_ts[features_ts_to_create], how = \"left\", left_index = True, right_index = True)\ndesc_continous_features(features_ts_to_create, train = t_)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T19:31:46.157181Z","iopub.execute_input":"2024-12-17T19:31:46.158117Z","iopub.status.idle":"2024-12-17T19:32:04.412356Z","shell.execute_reply.started":"2024-12-17T19:31:46.158078Z","shell.execute_reply":"2024-12-17T19:32:04.411297Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**q95_max_light_school** : the fewer children go to school and see daylight during school time, the more parents judge them as addicted to the internet !","metadata":{}},{"cell_type":"code","source":"t_ = train.merge(train_ts[features_ts_to_create], how = \"left\", left_index = True, right_index = True)\nother_plot = []\nif \"anglez_evening_kurt\" in features_ts_to_create:\n    t_[\"anglez_evening_kurt_clipped\"] = t_[\"anglez_evening_kurt\"].clip(0, 1)\n    other_plot.append(\"anglez_evening_kurt_clipped\")\nif \"q90_max_light_wevening\" in features_ts_to_create:\n    t_[\"q90_max_light_wevening_clipped\"] = t_[\"q90_max_light_wevening\"].clip(0, 400)\n    other_plot.append(\"q90_max_light_wevening_clipped\")\nif \"norm_yz_night_q75\" in features_ts_to_create:\n    t_[\"norm_yz_night_q75_clipped\"] = t_[\"norm_yz_night_q75\"].clip(0.8, 2)\n    other_plot.append(\"norm_yz_night_q75_clipped\")\nif \"norm_xz_evening_std\" in features_ts_to_create:\n    t_[\"norm_xz_evening_std_clipped\"] = t_[\"norm_xz_evening_std\"].clip(0.05, 0.5)\n    other_plot.append(\"norm_xz_evening_std_clipped\")\ndesc_continous_features(other_plot, train = t_)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T19:32:04.414159Z","iopub.execute_input":"2024-12-17T19:32:04.414468Z","iopub.status.idle":"2024-12-17T19:32:09.763601Z","shell.execute_reply.started":"2024-12-17T19:32:04.414437Z","shell.execute_reply":"2024-12-17T19:32:09.762450Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Debugging my test predictions : my public score is so low...\n# No good line below, except for debugging and to ensure predictions on test are made by the same ways than on oof samples\n'''\noof_scores = []; pred_scores = []\ndf = pd.DataFrame(index=train.dropna(subset=target).index)\nfor i in range(n_repeats):\n    a = pd.read_csv(f\"{output_path}/oof_label_xgb24_r{i+1}.csv\", index_col = \"id\")\n    df[i] = pd.read_csv(f\"{output_path}/pred_label_xgb24_r{i+1}.csv\", index_col = \"id\")\n    oof_scores.append(score_(train.dropna(subset=target)[target], a))\n    pred_scores.append(score_(train.dropna(subset=target)[target], df[i]))\n    print(f\"{oof_scores[-1:][0]:.3f} | {pred_scores[-1:][0]:.3f}\")\n\nprint(f\"{np.mean(oof_scores):.3f}(±{np.std(oof_scores):.3f}) {np.mean(pred_scores):.3f}(±{np.std(pred_scores):.3f}) \")\ntarget_values = sorted(list(train[target].unique()))\ncols = list(df.columns)\nfor i in target_values:\n    df[f\"nb_{i}\"] = 0\n    for s in cols:\n        df.loc[df[s] == i, f\"nb_{i}\"] += 1\n            \ndf[target] = np.argmax(df[[f\"nb_{i}\" for i in target_values]], axis = 1)\nprint(f\"{score_(train.dropna(subset=target)[target], df[target]):.3f}\")\n\n(df[9]!=a[target]).sum()\n'''","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T20:04:25.535068Z","iopub.execute_input":"2024-12-17T20:04:25.535459Z","iopub.status.idle":"2024-12-17T20:04:25.543009Z","shell.execute_reply.started":"2024-12-17T20:04:25.535426Z","shell.execute_reply":"2024-12-17T20:04:25.541959Z"},"_kg_hide-input":true,"_kg_hide-output":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Ensemble (Hard Vote)\nTraining an ensemble using majority vote will reduce the variance of predictions and the variance of the score.\n## Usefull functions","metadata":{}},{"cell_type":"code","source":"def score_after_hard_vote(mods_, df):\n\n    vs = sorted(list(df[target].unique()))\n    \n    for i in vs:\n        df[f\"nb_{i}\"] = 0\n        for s in mods_[1:][0]:\n            df.loc[df[s] == i, f\"nb_{i}\"] += 1\n            \n    df[\"res\"] = np.argmax(df[[f\"nb_{i}\" for i in vs]], axis = 1)\n    \n    return mods_[0], score_(df[target], df[\"res\"])\n    \n    \ndef fun_combinations(length, n_comb = n_repeats, replace = False, random_state = seed):\n    \n    if n_comb > comb(n_repeats, length):\n        print(f\"Warning : will return only {comb(n_repeats, length)} comb\")\n        n_comb = comb(n_repeats, length)\n    combs = []\n    np.random.seed(random_state)\n    while len(combs) < n_comb:\n        experience = sorted(np.random.choice(list(range(n_repeats)), size = length, replace = replace).tolist())\n        if experience not in combs:\n            combs.append(experience)\n    return combs\nprint(fun_combinations(3, n_comb = 10))\nprint(comb(n_repeats, 3))\nprint(factorial(n_repeats)/(factorial(3)* factorial(n_repeats-3)))\n\ndef run_experiences(experiences, title, do_plot = False):\n    \n    my_func = partial(score_after_hard_vote, df = df_oof_label)\n    all_scores = df_parallelize_run(my_func, experiences)\n    \n    sc = {}\n    for s in all_scores:\n        if s[0] not in list(sc.keys()):\n            sc[s[0]] = []\n        sc[s[0]].append(s[1])\n    \n    res = pd.DataFrame(sc).T\n    cols = list(res.columns)\n    res[\"mean\"] = res[cols].mean(axis=1)\n    res[\"std\"] = res[cols].std(axis=1)\n    res[\"rank1\"] = res[\"mean\"].rank(ascending = False)\n\n    ascending, higher_ = True, \"higher\"\n\n    _col_order = list(res[\"mean\"].T.sort_values(ascending = ascending).index)\n    res = res.loc[_col_order]\n    \n    if do_plot:\n        _t = list(res[\"mean\"].values)\n        _s = list(res[\"std\"].values)\n        _labels = [f\"{l:15} {v:.4f} (±{s:.4f})\" for l, v, s in zip(_col_order, _t, _s)]\n\n        fig, ax = plt.subplots(1, 1, figsize = (10, max(1, int(res.shape[0] * 1/2))))\n        ax.boxplot(res.T.loc[cols], vert = False)\n        ax.set_yticklabels(_labels)\n        ax.set_title(title)\n        ax.set_xlabel(f\"{metric_label} ({higher_} is better)\")\n        ax.xaxis.set_ticks_position(\"top\")\n        ax.xaxis.set_label_position('top')\n        ax.spines[[\"right\", \"bottom\"]].set_visible(False)\n    \n    return res[[\"mean\", \"std\"]]","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T20:05:22.474208Z","iopub.execute_input":"2024-12-17T20:05:22.475261Z","iopub.status.idle":"2024-12-17T20:05:22.492036Z","shell.execute_reply.started":"2024-12-17T20:05:22.475219Z","shell.execute_reply":"2024-12-17T20:05:22.490919Z"},"_kg_hide-input":true,"_kg_hide-output":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## Get oof data\nLet's get all prediction labels, from oof files and test prediction files :","metadata":{}},{"cell_type":"code","source":"train = train.dropna(subset = target)\n\ndf_oof_label, df_pred_label = pd.DataFrame(index=train.index), pd.DataFrame(index=test.index)\nfor r in range(n_repeats):\n    \n    for i, m in enumerate(list(all_results.keys())):\n        df_oof_label = pd.concat([\n            df_oof_label,\n            pd.read_csv(f\"{output_path}/oof_label_{m}_r{r+1}.csv\", index_col = \"id\").rename(columns={target:f\"{m}_{r}\"})\n        ], axis = 1)\n        \n        df_pred_label = pd.concat([\n            df_pred_label,\n            pd.read_csv(f\"{output_path}/pred_label_{m}_r{r+1}.csv\", index_col = \"id\").rename(columns={target:f\"{m}_{r}\"})\n        ], axis = 1)\n        \ndf_oof_label[target] = train[target]\ndisplay(df_oof_label.head())\ndisplay(df_pred_label.head())\n\nif not keep_files:\n    for root, dirs, files in os.walk(f\"./{output_path}/\", topdown=False):\n        for name in files:\n            os.remove(os.path.join(root, name))\n    os.rmdir(os.path.join(\"./\", f\"{output_path}\"))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T20:05:38.918371Z","iopub.execute_input":"2024-12-17T20:05:38.919444Z","iopub.status.idle":"2024-12-17T20:05:39.041765Z","shell.execute_reply.started":"2024-12-17T20:05:38.919387Z","shell.execute_reply":"2024-12-17T20:05:39.040573Z"},"_kg_hide-input":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## Best single model \nWhich model is the best without any stacking","metadata":{}},{"cell_type":"code","source":"models = list(all_params.keys())\nn_comb = 40 if n_repeats == 100 else 10","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T20:06:12.883278Z","iopub.execute_input":"2024-12-17T20:06:12.885078Z","iopub.status.idle":"2024-12-17T20:06:12.892350Z","shell.execute_reply.started":"2024-12-17T20:06:12.885004Z","shell.execute_reply":"2024-12-17T20:06:12.890879Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"all_experiences = []\nexperiences = fun_combinations(1, replace = False, random_state = seed, n_comb = n_comb)\nfor m in models:\n    for e in experiences:\n        all_experiences.append([f\"1 * {m}\", [f\"{m}_{e}\" for e in e]])\n            \nres = run_experiences(all_experiences, f\"{metric_label} of individual oof files\", do_plot = True)\nres","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T20:06:15.281749Z","iopub.execute_input":"2024-12-17T20:06:15.282135Z","iopub.status.idle":"2024-12-17T20:06:15.912966Z","shell.execute_reply.started":"2024-12-17T20:06:15.282100Z","shell.execute_reply":"2024-12-17T20:06:15.911376Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## Hard vote\n* Let's take randomly 1 *xgb24* oof files, compute QWK score, do it 40 times, compute the mean score and its standard deviation.\n* Let's take randomly 3 *xgb24* oof files, let them vote, compute QWK score, do it 40 times, compute the mean score and its standard deviation.\n* Let's take randomly 5 *xgb24* oof files, let them vote, compute QWK score, do it 40 times, compute the mean score and its standard deviation.\n* ...\n* Let's take randomly 97 *xgb24* oof files, let them vote, compute QWK score, do it 40 times, compute the mean score and its standard deviation.\n* and do it for the others GBM\n\n\nStandard deviation below is **underestimated** on the right side : the same prediction appears more frequently in the mix, thus the standard deviation decreases excessively. But the average increases !","metadata":{}},{"cell_type":"code","source":"all_experiences = []\nfor i in range(1, n_repeats - 2, 2):\n    experiences = fun_combinations(i, replace = False, random_state = seed, n_comb = n_comb)\n    for m in models:\n        for e in experiences:\n            all_experiences.append([f\"{i} * {m}\", [f\"{m}_{e_}\" for e_ in e]])\n            \nres = run_experiences(all_experiences, f\"{metric_label} of hard votes from 1 to {n_repeats-2} oof files\")\n\nres = res[[\"mean\", \"std\"]]\nt_ = pd.DataFrame(\n    res.index.to_series().str.split().to_list(), \n    columns = [\"n_repeats\", \"_\", \"gbm\"], index = res.index)[[\"n_repeats\", \"gbm\"]]\nt_[\"n_repeats\"] = t_[\"n_repeats\"].astype(np.int16)\nres = pd.concat([res, t_], axis=1)\n\nfig, ax = plt.subplots(len(models), 1, figsize = (15, 5 * len(models)))\nif len(models)==1: ax = [ax]\n\nfor i, m in enumerate(models):\n    df_toplot = res.loc[res[\"gbm\"]== m, [\"mean\", \"std\", \"n_repeats\"]].reset_index(drop=True).set_index(\"n_repeats\").sort_index()\n\n    q05 = df_toplot[\"mean\"] - 1.96 * df_toplot[\"std\"]\n    q95 = df_toplot[\"mean\"] + 1.96 * df_toplot[\"std\"]\n    \n    ax[i].plot(df_toplot[\"mean\"], linewidth=2) #mean curve.\n    ax[i].fill_between(df_toplot.index, q05, q95, color='b', alpha=.1)\n    ax[i].set_xlabel(f\"N distinct {m} during hard vote\")\n    ax[i].set_ylabel(\"QWK\")\n    ax[i].set_title(f\"QWK mean and confidence interval for each {n_comb} experiences of hard votes of {m}\")\n    ax[i].xaxis.set_ticks_position(\"top\")\n    ax[i].xaxis.set_label_position('top')\n    ax[i].spines[[\"right\", \"bottom\"]].set_visible(False)\n    \n    ax[i].axline(\n        (df_toplot.index.min(), df_toplot[\"mean\"].max()), \n        (df_toplot.index.max(), df_toplot[\"mean\"].max()),\n        color=\"g\",\n    )\n    mask = df_toplot[\"mean\"] == df_toplot[\"mean\"].max()\n    ind, mu, sigma = df_toplot.loc[mask].index.values[0], df_toplot.loc[mask, \"mean\"].values[0], df_toplot.loc[mask, \"std\"].values[0]\n    ax[i].plot(df_toplot.loc[mask].index.values[0], df_toplot[\"mean\"].max(), color=\"g\", marker=\"o\")\n    ax[i].text(df_toplot.loc[mask].index.values[0]-5, df_toplot[\"mean\"].max()+.001, f\"Max mean ({ind}) {mu:.4f} ±{sigma:.4f}\")\n    \n    mask = df_toplot[\"std\"] == df_toplot[\"std\"].min()\n    ind, mu, sigma = df_toplot.loc[mask].index.values[0], df_toplot.loc[mask, \"mean\"].values[0], df_toplot.loc[mask, \"std\"].values[0]\n    ax[i].plot(df_toplot.loc[mask].index.values[0], df_toplot.loc[mask, \"mean\"], color=\"b\", marker=\"o\")\n    ax[i].text(df_toplot.loc[mask].index.values[0]-5, df_toplot.loc[mask, \"mean\"]-.002, f\"Min std ({ind}) {mu:.4f} ±{sigma:.4f}\")\n    \n    mask = df_toplot.index == 1\n    ind, mu, sigma = 1, df_toplot.loc[mask, \"mean\"].values[0], df_toplot.loc[mask, \"std\"].values[0]\n    ax[i].plot(df_toplot.loc[mask].index.values[0], df_toplot.loc[mask, \"mean\"], color=\"r\", marker=\"o\")\n    ax[i].text(df_toplot.loc[mask].index.values[0], df_toplot.loc[mask, \"mean\"]-.002, f\"No hard vote ({ind}): {mu:.4f} ±{sigma:.4f}\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T20:06:24.733499Z","iopub.execute_input":"2024-12-17T20:06:24.733947Z","iopub.status.idle":"2024-12-17T20:06:25.316090Z","shell.execute_reply.started":"2024-12-17T20:06:24.733911Z","shell.execute_reply":"2024-12-17T20:06:25.314820Z"},"_kg_hide-input":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"best_result = res.loc[res[\"mean\"]==res[\"mean\"].max()].reset_index(drop=True).set_index(\"gbm\").to_dict(orient='index')\nmod_n1 = list(best_result.keys())[0]\nprint(best_result)\n\ndf_toplot = res.loc[res[\"gbm\"]== mod_n1, [\"mean\", \"std\", \"n_repeats\"]].reset_index(drop=True).set_index(\"n_repeats\").sort_index()\n\nfig, ax = plt.subplots(1, 1, figsize = (10, 5))\nx = np.linspace(.47, .5, 100)\n\nmask = df_toplot[\"std\"] == df_toplot[\"std\"].min()\n#print(df_toplot.loc[mask].index)\nmu, sigma = df_toplot.loc[mask, \"mean\"].values[0], df_toplot.loc[mask, \"std\"].values[0]\nax.plot(x, stats.norm.pdf(x, mu, sigma), label = f\"Lowest std QWK ({df_toplot.loc[mask].index[0]} repet.)\")\n\nmask = df_toplot[\"mean\"] == df_toplot[\"mean\"].max()\n#print(df_toplot.loc[mask].index)\nmu, sigma = df_toplot.loc[mask, \"mean\"].values[0], df_toplot.loc[mask, \"std\"].values[0]\nax.plot(x, stats.norm.pdf(x, mu, sigma), label = f\"Highest mean QWK ({df_toplot.loc[mask].index[0]} repet.)\")\n\nmask = df_toplot[\"mean\"] == df_toplot[\"mean\"].max()\n#print(df_toplot.loc[mask].index)\nmu, sigma = df_toplot.loc[mask, \"mean\"].values[0], df_toplot.loc[mask, \"std\"].values[0] * 3\nax.plot(x, stats.norm.pdf(x, mu, sigma), label = f\"Highest mean QWK ({df_toplot.loc[mask].index[0]} repet.) - stddev X 3\")\n\nmask = df_toplot.index == 1\n#print(df_toplot.loc[mask].index)\nmu, sigma = df_toplot.loc[mask, \"mean\"].values[0], df_toplot.loc[mask, \"std\"].values[0]\nax.plot(x, stats.norm.pdf(x, mu, sigma), label = f\"QWK without upvote ({df_toplot.loc[mask].index[0]} repet.)\")\n\nax.set_xlabel(\"QWK\")\nax.set_title(f\"Distribution of QWK with hard upvotes of {mod_n1}\")\nax.spines[[\"right\", \"top\"]].set_visible(False)\nax.legend();","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T20:14:35.190465Z","iopub.execute_input":"2024-12-17T20:14:35.190950Z","iopub.status.idle":"2024-12-17T20:14:35.470188Z","shell.execute_reply.started":"2024-12-17T20:14:35.190905Z","shell.execute_reply":"2024-12-17T20:14:35.469083Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"Take care : Std dev is underestimated below.\n# Final inference & submission","metadata":{}},{"cell_type":"code","source":"sub = pd.DataFrame(index = test.index)\ncomb_for_sub = fun_combinations(best_result[mod_n1][\"n_repeats\"], replace = False, random_state = seed, n_comb = 1)\nsub = pd.concat([sub, df_pred_label[[f'{mod_n1}_{v}' for v in comb_for_sub[0]]]], axis=1)\n\ntarget_values = sorted(list(train[target].unique()))\ncols = list(sub.columns)\nfor i in target_values:\n    sub[f\"nb_{i}\"] = 0\n    for s in cols:\n        sub.loc[sub[s] == i, f\"nb_{i}\"] += 1\n            \nsub[target] = np.argmax(sub[[f\"nb_{i}\" for i in target_values]], axis = 1)\nsub[target].to_csv(\"submission.csv\")\ndisplay(pd.read_csv(\"submission.csv\").head(20))\ndisplay(sub.head())\n# Not a good line below, except for debugging and to ensure predictions on test are made by the same ways than on oof samples\n#print(f\"Debugging test predictions : {score_(sub[target], train.dropna(subset=target)[target]):.3f}\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-17T20:14:01.834560Z","iopub.execute_input":"2024-12-17T20:14:01.835521Z","iopub.status.idle":"2024-12-17T20:14:01.902131Z","shell.execute_reply.started":"2024-12-17T20:14:01.835464Z","shell.execute_reply":"2024-12-17T20:14:01.901106Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}