{"metadata":{"kernelspec":{"name":"python3","display_name":"Python 3","language":"python"},"language_info":{"name":"python","version":"3.10.12","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":30786,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport os\nfrom sklearn.base import clone\nfrom sklearn.metrics import cohen_kappa_score, make_scorer, confusion_matrix\nfrom sklearn.model_selection import StratifiedKFold, KFold, train_test_split\nfrom sklearn.preprocessing import StandardScaler\nfrom sklearn.decomposition import PCA\nfrom scipy.optimize import minimize\nfrom scipy import stats\nfrom concurrent.futures import ThreadPoolExecutor\nfrom tqdm import tqdm\nimport warnings\nfrom sklearn.linear_model import ElasticNetCV, LassoCV, Lasso, LinearRegression\nfrom sklearn.ensemble import ExtraTreesRegressor, RandomForestRegressor\nfrom lightgbm import LGBMRegressor\nfrom xgboost import XGBRegressor\nfrom catboost import CatBoostRegressor\nimport optuna\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport numpy as np\nimport random\n\nwarnings.filterwarnings('ignore')","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:27:33.484826Z","iopub.execute_input":"2025-01-02T10:27:33.485178Z","iopub.status.idle":"2025-01-02T10:27:33.492558Z","shell.execute_reply.started":"2025-01-02T10:27:33.485152Z","shell.execute_reply":"2025-01-02T10:27:33.491326Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"SEED = 643\nn_splits = 10\noptimize_params = False\nn_trials = 25 # n_trials for optuna \nvoting = True\nbase_thresholds = [30, 50, 80]","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:27:33.506230Z","iopub.execute_input":"2025-01-02T10:27:33.506605Z","iopub.status.idle":"2025-01-02T10:27:33.520164Z","shell.execute_reply.started":"2025-01-02T10:27:33.506574Z","shell.execute_reply":"2025-01-02T10:27:33.518756Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"TRAIN_PATH = '/kaggle/input/child-mind-institute-problematic-internet-use/train.csv'\nTEST_PATH = '/kaggle/input/child-mind-institute-problematic-internet-use/test.csv'\nTRAIN_TS_PATH = '/kaggle/input/child-mind-institute-problematic-internet-use/series_train.parquet'\nTEST_TS_PATH = '/kaggle/input/child-mind-institute-problematic-internet-use/series_test.parquet'\nSUBMISSION_PATH = '/kaggle/input/child-mind-institute-problematic-internet-use/sample_submission.csv'\nOUTPUT_PATH = '/kaggle/working/'","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-01-02T10:27:33.521698Z","iopub.execute_input":"2025-01-02T10:27:33.522055Z","iopub.status.idle":"2025-01-02T10:27:33.541607Z","shell.execute_reply.started":"2025-01-02T10:27:33.522015Z","shell.execute_reply":"2025-01-02T10:27:33.540267Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Load datasets\ntrain = pd.read_csv(TRAIN_PATH)\ntest = pd.read_csv(TEST_PATH)","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:27:33.543556Z","iopub.execute_input":"2025-01-02T10:27:33.544014Z","iopub.status.idle":"2025-01-02T10:27:33.608513Z","shell.execute_reply.started":"2025-01-02T10:27:33.543970Z","shell.execute_reply":"2025-01-02T10:27:33.607292Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def perform_pca(train, test, n_components=None, random_state=42):\n    \"\"\"\n    Performs PCA on train, describes the Explained Variance Ratio and transforms train and test\n    \n    Returns: transformed train, transformed test, pca object\n    \n    \"\"\"\n    pca = PCA(n_components=n_components, random_state=random_state)\n    train_pca = pca.fit_transform(train)\n    test_pca = pca.transform(test)\n    \n    explained_variance_ratio = pca.explained_variance_ratio_\n    print(f\"Explained variance ratio of the components:\\n {explained_variance_ratio}\")\n    print(np.sum(explained_variance_ratio))\n    \n    train_pca_df = pd.DataFrame(train_pca, columns=[f'PC_{i+1}' for i in range(train_pca.shape[1])])\n    test_pca_df = pd.DataFrame(test_pca, columns=[f'PC_{i+1}' for i in range(test_pca.shape[1])])\n    \n    return train_pca_df, test_pca_df, pca","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:27:33.610042Z","iopub.execute_input":"2025-01-02T10:27:33.610354Z","iopub.status.idle":"2025-01-02T10:27:33.617695Z","shell.execute_reply.started":"2025-01-02T10:27:33.610328Z","shell.execute_reply":"2025-01-02T10:27:33.616251Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def time_features(df):\n    \"\"\"Function extracting Features from ActiGraph data of an individual\"\"\"\n    # Convert time_of_day to hours\n    df[\"hours\"] = df[\"time_of_day\"] // (3_600 * 1_000_000_000)\n    # Basic features \n    features = [\n        df[\"non-wear_flag\"].mean(),\n        df[\"enmo\"][df[\"enmo\"] >= 0.05].sum(),\n    ]\n    \n    # Define conditions for night, day, and no mask (full data)\n    night = ((df[\"hours\"] >= 22) | (df[\"hours\"] <= 5))\n    day = ((df[\"hours\"] <= 20) & (df[\"hours\"] >= 7))\n    no_mask = np.ones(len(df), dtype=bool)\n    \n    # List of columns of interest and masks\n    keys = [\"enmo\", \"anglez\", \"light\", \"battery_voltage\"]\n    masks = [no_mask, night, day]\n    \n    # Helper function for feature extraction\n    def extract_stats(data):\n        return [\n            data.mean(), \n            data.std(), \n            data.max(), \n            data.min(), \n            data.diff().mean(), \n            data.diff().std()\n        ]\n    \n    # Iterate over keys and masks to generate the statistics\n    for key in keys:\n        for mask in masks:\n            filtered_data = df.loc[mask, key]\n            features.extend(extract_stats(filtered_data))\n\n    return features\n\n# Code for parallelized computation of time series data from: Sheikh Muhammad Abdullah \n# https://www.kaggle.com/code/abdmental01/cmi-best-single-model\ndef process_file(filename, dirname):\n    # Process file and extract time features\n    df = pd.read_parquet(os.path.join(dirname, filename, 'part-0.parquet'))\n    df.drop('step', axis=1, inplace=True)\n    return time_features(df), filename.split('=')[1]\n\ndef load_time_series(dirname) -> pd.DataFrame:\n    # Load time series from directory in parallel\n    ids = os.listdir(dirname)\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, indexes = zip(*results)\n    \n    df = pd.DataFrame(stats, columns=[f\"stat_{i}\" for i in range(len(stats[0]))])\n    df['id'] = indexes\n    \n    return df","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:27:33.618468Z","iopub.execute_input":"2025-01-02T10:27:33.618912Z","iopub.status.idle":"2025-01-02T10:27:33.642136Z","shell.execute_reply.started":"2025-01-02T10:27:33.618785Z","shell.execute_reply":"2025-01-02T10:27:33.641127Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def feature_engineering(df):\n    season_cols = [col for col in df.columns if 'Season' in col]\n    df = df.drop(season_cols, axis=1) \n    \n    # From here on own features\n    def assign_group(age):\n        thresholds = [5, 6, 7, 8, 10, 12, 14, 17, 22]\n        for i, j in enumerate(thresholds):\n            if age <= j:\n                return i\n        return np.nan\n    \n    # Age groups\n    df[\"group\"] = df['Basic_Demos-Age'].apply(assign_group)\n    \n    # BMI \n    BMI_map = {0: 16.3,1: 15.9,2: 16.1,3: 16.8,4: 17.3,5: 19.2,6: 20.2,7: 22.3, 8: 23.6}\n    df['BMI_mean_norm'] = df[['Physical-BMI', 'BIA-BIA_BMI']].mean(axis=1) / df[\"group\"].map(BMI_map)\n    \n    # FGC zone aggregate\n    zones = ['FGC-FGC_CU_Zone', 'FGC-FGC_GSND_Zone', 'FGC-FGC_GSD_Zone',\n             'FGC-FGC_PU_Zone', 'FGC-FGC_SRL_Zone', 'FGC-FGC_SRR_Zone',\n             'FGC-FGC_TL_Zone']\n    \n    df['FGC_Zones_mean'] = df[zones].mean(axis=1)\n    df['FGC_Zones_min'] = df[zones].min(axis=1)\n    df['FGC_Zones_max'] = df[zones].max(axis=1)\n    \n    # Grip\n    GSD_max_map = {0: 9, 1: 9, 2: 9, 3: 9, 4: 16.2, 5: 19.9, 6: 26.1, 7: 31.3, 8: 35.4}\n    GSD_min_map = {0: 9, 1: 9, 2: 9, 3: 9, 4: 14.4, 5: 17.8, 6: 23.4, 7: 27.8, 8: 31.1}\n    \n    df['GS_max'] = df[['FGC-FGC_GSND', 'FGC-FGC_GSD']].max(axis=1) / df[\"group\"].map(GSD_max_map)\n    df['GS_min'] = df[['FGC-FGC_GSND', 'FGC-FGC_GSD']].min(axis=1) / df[\"group\"].map(GSD_min_map)\n    \n    # Curl-ups, push-ups, trunk-lifts... normalized based on age-group\n    cu_map = {0: 1.0, 1: 3.0, 2: 5.0, 3: 7.0, 4: 10.0, 5: 14.0, 6: 20.0, 7: 20.0, 8: 20.0}\n    pu_map = {0: 1.0, 1: 2.0, 2: 3.0, 3: 4.0, 4: 5.0, 5: 7.0, 6: 8.0, 7: 10.0, 8: 14.0}\n    tl_map = {0: 8.0, 1: 8.0, 2: 8.0, 3: 9.0, 4: 9.0, 5: 10.0, 6: 10.0, 7: 10.0, 8: 10.0}\n    \n    df[\"CU_norm\"] = df['FGC-FGC_CU'] / df['group'].map(cu_map)\n    df[\"PU_norm\"] = df['FGC-FGC_PU'] / df['group'].map(pu_map)\n    df[\"TL_norm\"] = df['FGC-FGC_TL'] / df['group'].map(tl_map)\n    \n    # Reach \n    df[\"SR_min\"] = df[['FGC-FGC_SRL', 'FGC-FGC_SRR']].min(axis=1)\n    df[\"SR_max\"] = df[['FGC-FGC_SRL', 'FGC-FGC_SRR']].max(axis=1)\n\n    # BIA Features\n    # Energy Expenditure\n    bmr_map = {0: 934.0, 1: 941.0, 2: 999.0, 3: 1048.0, 4: 1283.0, 5: 1255.0, 6: 1481.0, 7: 1519.0, 8: 1650.0}\n    dee_map = {0: 1471.0, 1: 1508.0, 2: 1640.0, 3: 1735.0, 4: 2132.0, 5: 2121.0, 6: 2528.0, 7: 2566.0, 8: 2793.0}\n    df[\"BMR_norm\"] = df[\"BIA-BIA_BMR\"] / df[\"group\"].map(bmr_map)\n    df[\"DEE_norm\"] = df[\"BIA-BIA_DEE\"] / df[\"group\"].map(dee_map)\n    df[\"DEE_BMR\"] = df[\"BIA-BIA_DEE\"] - df[\"BIA-BIA_BMR\"]\n\n    # FMM\n    ffm_map = {0: 42.0, 1: 43.0, 2: 49.0, 3: 54.0, 4: 60.0, 5: 76.0, 6: 94.0, 7: 104.0, 8: 111.0}\n    df[\"FFM_norm\"] = df[\"BIA-BIA_FFM\"] / df[\"group\"].map(ffm_map)\n\n    # ECW ICW\n    df[\"ICW_ECW\"] = df[\"BIA-BIA_ECW\"] / df[\"BIA-BIA_ICW\"]\n    \n    drop_feats = ['FGC-FGC_GSND', 'FGC-FGC_GSD', 'FGC-FGC_CU_Zone', 'FGC-FGC_GSND_Zone', 'FGC-FGC_GSD_Zone',\n                  'FGC-FGC_PU_Zone', 'FGC-FGC_SRL_Zone', 'FGC-FGC_SRR_Zone', 'FGC-FGC_TL_Zone',\n                  'Physical-BMI', 'BIA-BIA_BMI', 'FGC-FGC_CU', 'FGC-FGC_PU', 'FGC-FGC_TL', 'FGC-FGC_SRL', 'FGC-FGC_SRR',\n                 'BIA-BIA_BMR', 'BIA-BIA_DEE', 'BIA-BIA_Frame_num', \"BIA-BIA_FFM\"]\n    df = df.drop(drop_feats, axis=1) \n    return df","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:27:33.643204Z","iopub.execute_input":"2025-01-02T10:27:33.643610Z","iopub.status.idle":"2025-01-02T10:27:33.666215Z","shell.execute_reply.started":"2025-01-02T10:27:33.643572Z","shell.execute_reply":"2025-01-02T10:27:33.664751Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"train_ts = load_time_series(TRAIN_TS_PATH)\ntest_ts = load_time_series(TEST_TS_PATH)\n\ndf_train = train_ts.drop('id', axis=1)\ndf_test = test_ts.drop('id', axis=1)\n\n# Scaling before PCA\nscaler = StandardScaler()\ndf_train = pd.DataFrame(scaler.fit_transform(df_train), columns=df_train.columns)\ndf_test = pd.DataFrame(scaler.transform(df_test), columns=df_test.columns)\n\n# Mean imputation of actigraph data before PCA\nfor c in df_train.columns:\n    m = np.mean(df_train[c])\n    df_train[c].fillna(m, inplace=True)\n    df_test[c].fillna(m, inplace=True)\n\nprint(df_train.shape)\n\ndf_train_pca, df_test_pca, pca = perform_pca(df_train, df_test, n_components=15, random_state=SEED)\n\ndf_train_pca['id'] = train_ts['id']\ndf_test_pca['id'] = test_ts['id']\n\ntrain = pd.merge(train, df_train_pca, how=\"left\", on='id')\ntest = pd.merge(test, df_test_pca, how=\"left\", on='id')\ntrain.shape","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:27:33.669443Z","iopub.execute_input":"2025-01-02T10:27:33.669923Z","iopub.status.idle":"2025-01-02T10:29:05.812886Z","shell.execute_reply.started":"2025-01-02T10:27:33.669866Z","shell.execute_reply":"2025-01-02T10:29:05.809275Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def clean_features(df):\n    # Remove highly implausible values\n\n    # Clip Grip\n    df[['FGC-FGC_GSND', 'FGC-FGC_GSD']] = df[['FGC-FGC_GSND', 'FGC-FGC_GSD']].clip(lower=9, upper=60)\n    # Remove implausible body-fat\n    df[\"BIA-BIA_Fat\"] = np.where(df[\"BIA-BIA_Fat\"] < 5, np.nan, df[\"BIA-BIA_Fat\"])\n    df[\"BIA-BIA_Fat\"] = np.where(df[\"BIA-BIA_Fat\"] > 60, np.nan, df[\"BIA-BIA_Fat\"])\n    # Basal Metabolic Rate\n    df[\"BIA-BIA_BMR\"] = np.where(df[\"BIA-BIA_BMR\"] > 4000, np.nan, df[\"BIA-BIA_BMR\"])\n    # Daily Energy Expenditure\n    df[\"BIA-BIA_DEE\"] = np.where(df[\"BIA-BIA_DEE\"] > 8000, np.nan, df[\"BIA-BIA_DEE\"])\n    # Bone Mineral Content\n    df[\"BIA-BIA_BMC\"] = np.where(df[\"BIA-BIA_BMC\"] <= 0, np.nan, df[\"BIA-BIA_BMC\"])\n    df[\"BIA-BIA_BMC\"] = np.where(df[\"BIA-BIA_BMC\"] > 10, np.nan, df[\"BIA-BIA_BMC\"])\n    # Fat Free Mass Index\n    df[\"BIA-BIA_FFM\"] = np.where(df[\"BIA-BIA_FFM\"] <= 0, np.nan, df[\"BIA-BIA_FFM\"])\n    df[\"BIA-BIA_FFM\"] = np.where(df[\"BIA-BIA_FFM\"] > 300, np.nan, df[\"BIA-BIA_FFM\"])\n    # Fat Mass Index\n    df[\"BIA-BIA_FMI\"] = np.where(df[\"BIA-BIA_FMI\"] < 0, np.nan, df[\"BIA-BIA_FMI\"])\n    # Extra Cellular Water\n    df[\"BIA-BIA_ECW\"] = np.where(df[\"BIA-BIA_ECW\"] > 100, np.nan, df[\"BIA-BIA_ECW\"])\n    # Intra Cellular Water\n    # df[\"BIA-BIA_ICW\"] = np.where(df[\"BIA-BIA_ICW\"] > 100, np.nan, df[\"BIA-BIA_ICW\"])\n    # Lean Dry Mass\n    df[\"BIA-BIA_LDM\"] = np.where(df[\"BIA-BIA_LDM\"] > 100, np.nan, df[\"BIA-BIA_LDM\"])\n    # Lean Soft Tissue\n    df[\"BIA-BIA_LST\"] = np.where(df[\"BIA-BIA_LST\"] > 300, np.nan, df[\"BIA-BIA_LST\"])\n    # Skeletal Muscle Mass\n    df[\"BIA-BIA_SMM\"] = np.where(df[\"BIA-BIA_SMM\"] > 300, np.nan, df[\"BIA-BIA_SMM\"])\n    # Total Body Water\n    df[\"BIA-BIA_TBW\"] = np.where(df[\"BIA-BIA_TBW\"] > 300, np.nan, df[\"BIA-BIA_TBW\"])\n    \n    return df\n\ntrain = clean_features(train)\ntest = clean_features(test)","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:05.814792Z","iopub.execute_input":"2025-01-02T10:29:05.815270Z","iopub.status.idle":"2025-01-02T10:29:05.875220Z","shell.execute_reply.started":"2025-01-02T10:29:05.815233Z","shell.execute_reply":"2025-01-02T10:29:05.873929Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# def assign_group(age):\n#     thresholds = [5, 6, 7, 8, 10, 12, 14, 17, 22]\n#     for i, j in enumerate(thresholds):\n#         if age <= j:\n#             return i\n#     return np.nan\n\n# fit = ['BIA-BIA_SMM',\n#        'BIA-BIA_TBW']\n\n# # train[['FGC-FGC_GSND', 'FGC-FGC_GSD']] = train[['FGC-FGC_GSND', 'FGC-FGC_GSD']].clip(lower=9, upper=60)\n# # train['GSD_max'] = train[['FGC-FGC_GSND', 'FGC-FGC_GSD']].max(axis=1)\n# # train['GSD_min'] = train[['FGC-FGC_GSND', 'FGC-FGC_GSD']].min(axis=1)\n# train[\"DEE_BMR\"] = train['BIA-BIA_DEE'] - train['BIA-BIA_BMR']\n# train[\"group\"] = train['Basic_Demos-Age'].apply(assign_group)\n# train.value_counts(\"group\")\n# np.round(train.groupby([\"group\"])[fit].mean())","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:05.876194Z","iopub.execute_input":"2025-01-02T10:29:05.876505Z","iopub.status.idle":"2025-01-02T10:29:05.881061Z","shell.execute_reply.started":"2025-01-02T10:29:05.876479Z","shell.execute_reply":"2025-01-02T10:29:05.879857Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# train[fit].describe()","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:05.882295Z","iopub.execute_input":"2025-01-02T10:29:05.882723Z","iopub.status.idle":"2025-01-02T10:29:05.907727Z","shell.execute_reply.started":"2025-01-02T10:29:05.882683Z","shell.execute_reply":"2025-01-02T10:29:05.906530Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# info = pd.read_csv(\"/kaggle/input/child-mind-institute-problematic-internet-use/data_dictionary.csv\")\n# info[33:]","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:05.908791Z","iopub.execute_input":"2025-01-02T10:29:05.909385Z","iopub.status.idle":"2025-01-02T10:29:05.931233Z","shell.execute_reply.started":"2025-01-02T10:29:05.909344Z","shell.execute_reply":"2025-01-02T10:29:05.929956Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Perform feature_engineering on train and test\ntrain = feature_engineering(train)\ntest = feature_engineering(test)","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:05.932482Z","iopub.execute_input":"2025-01-02T10:29:05.933145Z","iopub.status.idle":"2025-01-02T10:29:06.015791Z","shell.execute_reply.started":"2025-01-02T10:29:05.933106Z","shell.execute_reply":"2025-01-02T10:29:06.014662Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def bin_data(train, test, columns, n_bins=10):\n    # Combine train and test for consistent bin edges\n    combined = pd.concat([train, test], axis=0)\n    \n    bin_edges = {}\n    for col in columns:\n        # Compute quantile bin edges\n        edges = pd.qcut(combined[col], n_bins, retbins=True, labels=range(n_bins), duplicates=\"drop\")[1]\n        bin_edges[col] = edges\n    \n    # Apply the same bin edges to both train and test\n    for col, edges in bin_edges.items():\n        train[col] = pd.cut(\n            train[col], bins=edges, labels=range(len(edges) - 1), include_lowest=True\n        ).astype(float)\n        test[col] = pd.cut(\n            test[col], bins=edges, labels=range(len(edges) - 1), include_lowest=True\n        ).astype(float)\n    \n    return train, test\n\ncolumns_to_bin = [\n    \"PAQ_A-PAQ_A_Total\", \"BMR_norm\", \"DEE_norm\", \"GS_min\", \"GS_max\", \"BIA-BIA_FFMI\", \n    \"BIA-BIA_BMC\", \"Physical-HeartRate\", \"BIA-BIA_ICW\", \"Fitness_Endurance-Time_Sec\", \n    \"BIA-BIA_LDM\", \"BIA-BIA_SMM\", \"BIA-BIA_TBW\", \"DEE_BMR\", \"ICW_ECW\"\n]\n# Bin specified columns in train and test\ntrain, test = bin_data(train, test, columns_to_bin, n_bins=10)","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:06.016895Z","iopub.execute_input":"2025-01-02T10:29:06.017284Z","iopub.status.idle":"2025-01-02T10:29:06.089714Z","shell.execute_reply.started":"2025-01-02T10:29:06.017255Z","shell.execute_reply":"2025-01-02T10:29:06.088501Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Features to exclude, because they're not in test\nexclude = ['PCIAT-Season', 'PCIAT-PCIAT_01', 'PCIAT-PCIAT_02', 'PCIAT-PCIAT_03',\n           'PCIAT-PCIAT_04', 'PCIAT-PCIAT_05', 'PCIAT-PCIAT_06', 'PCIAT-PCIAT_07',\n           'PCIAT-PCIAT_08', 'PCIAT-PCIAT_09', 'PCIAT-PCIAT_10', 'PCIAT-PCIAT_11',\n           'PCIAT-PCIAT_12', 'PCIAT-PCIAT_13', 'PCIAT-PCIAT_14', 'PCIAT-PCIAT_15',\n           'PCIAT-PCIAT_16', 'PCIAT-PCIAT_17', 'PCIAT-PCIAT_18', 'PCIAT-PCIAT_19',\n           'PCIAT-PCIAT_20', 'PCIAT-PCIAT_Total', 'sii', 'id']\n\ny_model = \"PCIAT-PCIAT_Total\" # Score, target for the model\ny_comp = \"sii\" # Index, target of the competition\nfeatures = [f for f in train.columns if f not in exclude]\n\n# Categorical features\ncat_c = []\n\n# Mapping of categorical features (already dropped in feature engineering)\nfor col in cat_c:\n    a_map = {}\n    all_unique = set(train[col].unique()) | set(test[col].unique())\n    for i, value in enumerate(all_unique):\n        a_map[value] = i\n\n    train[col] = train[col].map(a_map)\n    test[col] = test[col].map(a_map)\n    \ntrain = train[train[\"sii\"].notna()] # Keep rows where target is available\ntrain.shape","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:06.091038Z","iopub.execute_input":"2025-01-02T10:29:06.091540Z","iopub.status.idle":"2025-01-02T10:29:06.110184Z","shell.execute_reply.started":"2025-01-02T10:29:06.091482Z","shell.execute_reply":"2025-01-02T10:29:06.107949Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Plot distribution of total scores which determine the sii\n# Note the excess zeros -> consider other objective functions\nsns.set_theme(style=\"whitegrid\")\nplt.hist(train['PCIAT-PCIAT_Total'], bins=50, color=\"darkorange\")\nplt.title('Score Distribution')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:06.111759Z","iopub.execute_input":"2025-01-02T10:29:06.112222Z","iopub.status.idle":"2025-01-02T10:29:06.576733Z","shell.execute_reply.started":"2025-01-02T10:29:06.112174Z","shell.execute_reply":"2025-01-02T10:29:06.575372Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"class Impute_With_Model:\n    \n    def __init__(self, na_frac=0.5, min_samples=0):\n        self.model_dict = {} # Dictionary storing models for imputation\n        self.mean_dict = {} # Dictionary storing mean of feature\n        self.features = None\n        self.na_frac = na_frac # Maximum fraction of missing values allowed for features to be used\n        self.min_samples = min_samples # Minimum number of samples required for model fitting\n        \n    def find_features(self, data, feature, tmp_features):\n        # Finds valid features where the fraction of missing values is <= na_frac\n        missing_rows = data[feature].isna()\n        na_fraction = data[missing_rows][tmp_features].isna().mean(axis=0)\n        valid_features = np.array(tmp_features)[na_fraction <= self.na_frac]\n        return valid_features\n\n    def fit_models(self, model, data, features):\n        self.features = features\n        n_data = data.shape[0]\n        for feature in features:\n            self.mean_dict[feature] = np.mean(data[feature])\n        # Iterate over all features\n        for feature in tqdm(features):\n            # Impute if there are missing values in the data\n            if data[feature].isna().sum() > 0:\n                model_clone = clone(model)\n                X = data[data[feature].notna()].copy() # Select data where target values are available as trainings data\n                tmp_features = [f for f in features if f != feature]\n                tmp_features = self.find_features(data, feature, tmp_features)\n                if len(tmp_features) >= 1 and X.shape[0] > self.min_samples:\n                    # Fit model if enough features and sufficient samples\n                    for f in tmp_features:\n                        X[f] = X[f].fillna(self.mean_dict[f])\n                    model_clone.fit(X[tmp_features], X[feature])\n                    # Add model and features to dictionary\n                    self.model_dict[feature] = (model_clone, tmp_features.copy())\n                else:\n                    # Revert to mean imputation if too few \n                    self.model_dict[feature] = (\"mean\", np.mean(data[feature]))\n            \n    def impute(self, data):\n        imputed_data = data.copy()\n        # Iterate over models\n        for feature, model in self.model_dict.items():\n            missing_rows = imputed_data[feature].isna() # Identify rows to be imputed\n            if missing_rows.any():\n                if model[0] == \"mean\":\n                    # Mean imputation if \"mean\"\n                    imputed_data[feature].fillna(model[1], inplace=True)\n                else:\n                    # Prepare data for imputation and predict\n                    tmp_features = [f for f in self.features if f != feature]\n                    X_missing = data.loc[missing_rows, tmp_features].copy()\n                    for f in tmp_features:\n                        X_missing[f] = X_missing[f].fillna(self.mean_dict[f])\n                    imputed_data.loc[missing_rows, feature] = model[0].predict(X_missing[model[1]])\n        return imputed_data","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:06.578500Z","iopub.execute_input":"2025-01-02T10:29:06.579000Z","iopub.status.idle":"2025-01-02T10:29:06.599666Z","shell.execute_reply.started":"2025-01-02T10:29:06.578950Z","shell.execute_reply":"2025-01-02T10:29:06.597302Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"missing = pd.DataFrame(train.isna().sum() / len(train))\nmissing[missing[0] > 0.3][:60]","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:06.604380Z","iopub.execute_input":"2025-01-02T10:29:06.604754Z","iopub.status.idle":"2025-01-02T10:29:06.647798Z","shell.execute_reply.started":"2025-01-02T10:29:06.604723Z","shell.execute_reply":"2025-01-02T10:29:06.646568Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"model = LassoCV(cv=5, random_state=SEED)\nimputer = Impute_With_Model(na_frac=0.4) \nimputer.fit_models(model, train, features)\ntrain = imputer.impute(train)\ntest = imputer.impute(test)","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:06.649760Z","iopub.execute_input":"2025-01-02T10:29:06.650197Z","iopub.status.idle":"2025-01-02T10:29:24.134499Z","shell.execute_reply.started":"2025-01-02T10:29:06.650147Z","shell.execute_reply":"2025-01-02T10:29:24.132597Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Code for finding optimal thresholds copied from: Michael Semenoff\n# https://www.kaggle.com/code/michaelsemenoff/cmi-actigraphy-feature-engineering-selection\ndef round_with_thresholds(raw_preds, thresholds):\n    return np.where(raw_preds < thresholds[0], int(0),\n                    np.where(raw_preds < thresholds[1], int(1),\n                             np.where(raw_preds < thresholds[2], int(2), int(3))))\n\ndef optimize_thresholds(y_true, raw_preds, start_vals=[0.5, 1.5, 2.5]):\n    def fun(thresholds, y_true, raw_preds):\n        rounded_preds = round_with_thresholds(raw_preds, thresholds)\n        return -cohen_kappa_score(y_true, rounded_preds, weights='quadratic')\n\n    res = minimize(fun, x0=start_vals, args=(y_true, raw_preds), method='Powell')\n    assert res.success\n    return res.x","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:24.136538Z","iopub.execute_input":"2025-01-02T10:29:24.137135Z","iopub.status.idle":"2025-01-02T10:29:24.146993Z","shell.execute_reply.started":"2025-01-02T10:29:24.137088Z","shell.execute_reply":"2025-01-02T10:29:24.145641Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def calculate_weights(series):\n    # Create bins for the target variable and assign weights based on frequency\n    bins = pd.cut(series, bins=10, labels=False)\n    weights = bins.value_counts().reset_index()\n    weights.columns = ['target_bins', 'count']\n    weights['count'] = 1 / weights['count']\n    weight_map = weights.set_index('target_bins')['count'].to_dict()\n    weights = bins.map(weight_map)\n    return weights / weights.mean() ","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:24.148251Z","iopub.execute_input":"2025-01-02T10:29:24.148690Z","iopub.status.idle":"2025-01-02T10:29:24.173316Z","shell.execute_reply.started":"2025-01-02T10:29:24.148650Z","shell.execute_reply":"2025-01-02T10:29:24.172147Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def cross_validate(model_, data, features, score_col, index_col, cv, sample_weights=False, verbose=False):\n    \"\"\"\n    Perform cross-validation with a given model and compute the out-of-fold \n    predictions and Cohen's Kappa score for each fold.\n\n    Returns:\n    float: Mean Kappa score across all folds.\n    array: Out-of-fold score predictions for the entire dataset.\n    \"\"\"\n    kappa_scores = [] \n    oof_score_predictions = np.zeros(len(data))  \n\n    score_to_index_thresholds = base_thresholds  \n    thresholds = []\n    for fold_idx, (train_idx, val_idx) in enumerate(cv.split(data, data[index_col])):\n        X_train, X_val = data[features].iloc[train_idx], data[features].iloc[val_idx]\n        y_train_score = data[score_col].iloc[train_idx] \n        y_train_index = data[index_col].iloc[train_idx]\n        y_val_score = data[score_col].iloc[val_idx]      \n        y_val_index = data[index_col].iloc[val_idx]     \n        \n        # Train model with sample weights if provided\n        if sample_weights:\n            weights = calculate_weights(y_train_score)\n            model_.fit(X_train, y_train_score, sample_weight=weights)\n        else:\n            model_.fit(X_train, y_train_score)\n\n        y_pred_train_score = model_.predict(X_train)\n        y_pred_val_score = model_.predict(X_val)\n        \n        oof_score_predictions[val_idx] = y_pred_val_score \n\n        # Find optimal threshold in sample \n        t_1 = optimize_thresholds(y_train_index, y_pred_train_score, start_vals=base_thresholds)\n        thresholds.append(t_1)\n\n        y_pred_val_index = round_with_thresholds(y_pred_val_score, t_1)\n\n        kappa_score = cohen_kappa_score(y_val_index, y_pred_val_index, weights='quadratic')\n        kappa_scores.append(kappa_score)\n        \n        if verbose:\n            print(f\"Fold {fold_idx}: Optimized Kappa Score = {kappa_score}\")\n    \n    if verbose:\n        print(f\"## Mean CV Kappa Score: {np.mean(kappa_scores)} ##\")\n        print(f\"## Std CV: {np.std(kappa_scores)}\")\n    \n    return np.mean(kappa_scores), oof_score_predictions, thresholds\n\ndef n_cross_validate(model_, data, features, score_col, index_col, cv, seeds, sample_weights=False, verbose=False):\n    # Performs repeated cross-validation by reseeding the cv object\n    scores = []\n    for seed in seeds:\n        cv.random_state=seed\n        score, oof, _ = cross_validate(model_, data, features, score_col, index_col, cv, sample_weights=True, verbose=False)\n        scores.append(score)\n    return score, oof","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:24.174611Z","iopub.execute_input":"2025-01-02T10:29:24.175039Z","iopub.status.idle":"2025-01-02T10:29:24.196421Z","shell.execute_reply.started":"2025-01-02T10:29:24.174997Z","shell.execute_reply":"2025-01-02T10:29:24.195157Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def objective(trial, model_type, X, features, score_col, index_col, cv, sample_weights=False):\n    # Parameter space to explore if model is xgboost\n    if model_type == 'xgboost':\n        params = {\n            'objective': trial.suggest_categorical('objective', ['reg:tweedie', 'reg:pseudohubererror']),\n            'random_state': SEED,\n            'num_parallel_tree': trial.suggest_int('num_parallel_tree', 2, 30),\n            'n_estimators': trial.suggest_int('n_estimators', 100, 300),\n            'max_depth': trial.suggest_int('max_depth', 2, 4),\n            'learning_rate': trial.suggest_loguniform('learning_rate', 0.02, 0.05),\n            'subsample': trial.suggest_float('subsample', 0.5, 0.8),\n            'colsample_bytree': trial.suggest_float('colsample_bytree', 0.5, 0.8),\n            'reg_alpha': trial.suggest_loguniform('reg_alpha', 1e-5, 1e-1),\n            'reg_lambda': trial.suggest_loguniform('reg_lambda', 1e-5, 1e-1),\n        }\n        if params['objective'] == 'reg:tweedie':\n            params['tweedie_variance_power'] = trial.suggest_float('tweedie_variance_power', 1, 2)\n        model = XGBRegressor(**params, use_label_encoder=False)\n    \n    # Parameter space to explore if model is lightgbm\n    elif model_type == 'lightgbm':\n        params = {\n            'objective': trial.suggest_categorical('objective', ['poisson', 'tweedie', 'regression']),\n            'random_state': SEED,\n            'verbosity': -1,\n            'n_estimators': trial.suggest_int('n_estimators', 100, 300),\n            'max_depth': trial.suggest_int('max_depth', 2, 4),\n            'learning_rate': trial.suggest_loguniform('learning_rate', 0.01, 0.05),\n            'subsample': trial.suggest_float('subsample', 0.5, 0.8),\n            'colsample_bytree': trial.suggest_float('colsample_bytree', 0.5, 0.8),\n            'min_data_in_leaf': trial.suggest_int('min_data_in_leaf', 20, 100)\n        }\n        if params['objective'] == 'tweedie':\n            params['tweedie_variance_power'] = trial.suggest_float('tweedie_variance_power', 1, 2)\n        model = LGBMRegressor(**params)\n    \n    # Parameter space to explore if model is catboost\n    elif model_type == 'catboost':\n        params = {\n            'loss_function': trial.suggest_categorical('objective', ['Tweedie:variance_power=1.5', \n                                                                     'Poisson', 'RMSE']),\n            'random_state': SEED,\n            'iterations': trial.suggest_int('iterations', 100, 300),\n            'depth': trial.suggest_int('depth', 2, 4),\n            'learning_rate': trial.suggest_loguniform('learning_rate', 0.01, 0.05),\n            'l2_leaf_reg': trial.suggest_loguniform('l2_leaf_reg', 1e-3, 1e-1),\n            'subsample': trial.suggest_float('subsample', 0.5, 0.7),\n            'bagging_temperature': trial.suggest_float('bagging_temperature', 0.0, 1.0),\n            'random_strength': trial.suggest_float('random_strength', 1e-3, 10.0),\n            'min_data_in_leaf': trial.suggest_int('min_data_in_leaf', 20, 60),\n        }\n        model = CatBoostRegressor(**params, verbose=0)\n    \n    else:\n        raise ValueError(f\"Unsupported model_type: {model_type}\")\n        \n    seeds = [random.randint(1, 10000) for _ in range(20)] # Seeds for repeated KFold\n\n    score, _ = n_cross_validate(model, X, features, score_col, index_col, cv, seeds, sample_weights=True, verbose=True)\n\n    return score\n\ndef run_optimization(X, features, score_col, index_col, model_type, n_trials=30, cv=None, sample_weights=False):\n    study = optuna.create_study(direction=\"maximize\")\n    study.optimize(lambda trial: objective(trial, model_type, X, features, score_col, index_col, cv, sample_weights), \n                   n_trials=n_trials)\n    \n    print(f\"Best params for {model_type}: {study.best_params}\")\n    print(f\"Best score: {study.best_value}\")\n    return study.best_params","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:24.198175Z","iopub.execute_input":"2025-01-02T10:29:24.198608Z","iopub.status.idle":"2025-01-02T10:29:24.226673Z","shell.execute_reply.started":"2025-01-02T10:29:24.198567Z","shell.execute_reply":"2025-01-02T10:29:24.225159Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Replace if subsets for features have been selected\n\n# List of manually selected features NOT to include\nexclude = [\n    \"PC_9\", \"PC_12\", \"Fitness_Endurance-Max_Stage\", \"Basic_Demos-Sex\", \"BMI_mean_norm\", \"PC_11\", \n    \"PC_8\", \"FGC_Zones_min\", \"Physical-Systolic_BP\", \"PC_4\", \"BIA-BIA_FMI\", \"BIA-BIA_LST\", \"Physical-Diastolic_BP\", \n    \"BIA-BIA_ECW\", \"Fitness_Endurance-Time_Mins\", \"PAQ_C-PAQ_C_Total\", \"PC_10\", \"BIA-BIA_Fat\", \"FFM_norm\", \"PC_14\", \"PC_7\"\n]\n\nreduced_features = [f for f in features if f not in exclude]\n\nlgb_features = reduced_features\nxgb_features = reduced_features\ncat_features = reduced_features\nprint(len(reduced_features))","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:24.228007Z","iopub.execute_input":"2025-01-02T10:29:24.228440Z","iopub.status.idle":"2025-01-02T10:29:24.254657Z","shell.execute_reply.started":"2025-01-02T10:29:24.228398Z","shell.execute_reply":"2025-01-02T10:29:24.253118Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Parameters for LGBM, XGB and CatBoost\nlgb_params = {\n    'objective': 'poisson', \n    'n_estimators': 295, \n    'max_depth': 4, \n    'learning_rate': 0.04505693066482616, \n    'subsample': 0.6042489155604022, \n    'colsample_bytree': 0.5021876720502726, \n    'min_data_in_leaf': 100\n}\n\nxgb_params = {\n    'objective': 'reg:tweedie', \n    'num_parallel_tree': 12, \n    'n_estimators': 236, \n    'max_depth': 3, \n    'learning_rate': 0.04223740904479563, \n    'subsample': 0.7157264603586825, \n    'colsample_bytree': 0.7897918901977528, \n    'reg_alpha': 0.005335705058190553, \n    'reg_lambda': 0.0001897435318347022, \n    'tweedie_variance_power': 1.1393958601390142\n}\n\nxgb_params_2 = {\n    'objective': 'reg:tweedie', \n    'num_parallel_tree': 18, \n    'n_estimators': 175, \n    'max_depth': 3, \n    'learning_rate': 0.032620453423049305, \n    'subsample': 0.6155579670568023, \n    'colsample_bytree': 0.5988773292417443, \n    'reg_alpha': 0.0028895066837627205, \n    'reg_lambda': 0.002232531512636924, \n    'tweedie_variance_power': 1.1708678482038286\n}\n\ncat_params = {\n    'objective': 'RMSE', \n    'iterations': 238, \n    'depth': 4, \n    'learning_rate': 0.044523361750173816, \n    'l2_leaf_reg': 0.09301285673435761, \n    'subsample': 0.6902492783438681, \n    'bagging_temperature': 0.3007304771330199, \n    'random_strength': 3.562201626987314, \n    'min_data_in_leaf': 60\n}\n\nxtrees_params = {\n    'n_estimators': 500, \n    'max_depth': 15, \n    'min_samples_leaf': 20, \n    'bootstrap': False\n}","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:24.255700Z","iopub.execute_input":"2025-01-02T10:29:24.255994Z","iopub.status.idle":"2025-01-02T10:29:24.278014Z","shell.execute_reply.started":"2025-01-02T10:29:24.255961Z","shell.execute_reply":"2025-01-02T10:29:24.276599Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"kf = StratifiedKFold(n_splits=n_splits, shuffle=True, random_state=SEED)","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:24.279209Z","iopub.execute_input":"2025-01-02T10:29:24.279602Z","iopub.status.idle":"2025-01-02T10:29:24.306107Z","shell.execute_reply.started":"2025-01-02T10:29:24.279560Z","shell.execute_reply":"2025-01-02T10:29:24.304904Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"if optimize_params:\n    # LightGBM Optimization\n    lgb_params = run_optimization(train, lgb_features, 'PCIAT-PCIAT_Total', 'sii', 'lightgbm', n_trials=n_trials, cv=kf, sample_weights=True)\n\n    # XGBoost Optimization\n    xgb_params = run_optimization(train, xgb_features, 'PCIAT-PCIAT_Total', 'sii', 'xgboost', n_trials=n_trials, cv=kf, sample_weights=True)\n\n    # CatBoost Optimization\n    cat_params = run_optimization(train, cat_features, 'PCIAT-PCIAT_Total', 'sii', 'catboost', n_trials=n_trials, cv=kf, sample_weights=True)","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:24.307318Z","iopub.execute_input":"2025-01-02T10:29:24.307726Z","iopub.status.idle":"2025-01-02T10:29:24.325936Z","shell.execute_reply.started":"2025-01-02T10:29:24.307684Z","shell.execute_reply":"2025-01-02T10:29:24.324593Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Define models\nlgb_model = LGBMRegressor(**lgb_params, random_state=SEED, verbosity=-1)\nxgb_model = XGBRegressor(**xgb_params, random_state=SEED, verbosity=0)\nxgb_model_2 = XGBRegressor(**xgb_params_2, random_state=SEED, verbosity=0)\ncat_model = CatBoostRegressor(**cat_params, random_state=SEED, verbose=0)\nxtrees_model = ExtraTreesRegressor(**xtrees_params, random_state=SEED)\n\nweights = calculate_weights(train['PCIAT-PCIAT_Total'])\n\n# Cross-validate LGBM model\nscore_lgb, oof_lgb, lgb_thresholds = cross_validate(\n    lgb_model, train, lgb_features, 'PCIAT-PCIAT_Total', 'sii', kf, verbose=True, sample_weights=True\n)\n# Fit final model and predict test samples\nlgb_model.fit(train[lgb_features], train['PCIAT-PCIAT_Total'], sample_weight=weights)\ntest_lgb = lgb_model.predict(test[lgb_features])\n\n# Cross-validate XGBoost model\nscore_xgb, oof_xgb, xgb_thresholds = cross_validate(\n    xgb_model, train, xgb_features, 'PCIAT-PCIAT_Total', 'sii', kf, verbose=True, sample_weights=True\n)\n# Fit final model and predict test samples\nxgb_model.fit(train[xgb_features], train['PCIAT-PCIAT_Total'], sample_weight=weights)\ntest_xgb = xgb_model.predict(test[xgb_features])\n\n# Cross-validate XGBoost model 2\nscore_xgb_2, oof_xgb_2, xgb_2_thresholds = cross_validate(\n    xgb_model_2, train, xgb_features, 'PCIAT-PCIAT_Total', 'sii', kf, verbose=True, sample_weights=True\n)\n# Fit final model and predict test samples\nxgb_model_2.fit(train[xgb_features], train['PCIAT-PCIAT_Total'], sample_weight=weights)\ntest_xgb_2 = xgb_model_2.predict(test[xgb_features])\n\n# Cross-validate CatBoost model\nscore_cat, oof_cat, cat_thresholds = cross_validate(\n    cat_model, train, cat_features, 'PCIAT-PCIAT_Total', 'sii', kf, verbose=True, sample_weights=True\n)\n# Fit final model and predict test samples\ncat_model.fit(train[cat_features], train['PCIAT-PCIAT_Total'], sample_weight=weights)\ntest_cat = cat_model.predict(test[cat_features])\n\n# Cross-validate ExtraTreesRegressor model\nscore_xtrees, oof_xtrees, xtrees_thresholds = cross_validate(\n    xtrees_model, train, reduced_features, 'PCIAT-PCIAT_Total', 'sii', kf, verbose=True, sample_weights=True\n)\n# Fit final model and predict test samples\nxtrees_model.fit(train[reduced_features], train['PCIAT-PCIAT_Total'], sample_weight=weights)\ntest_xtrees = xtrees_model.predict(test[reduced_features])\n\n# Print overall mean Kappa score for all models\nprint(f'Overall Mean Kappa: {np.mean([score_lgb, score_xgb, score_cat, score_xtrees])}') # Ensemble score likely higher","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:29:24.327433Z","iopub.execute_input":"2025-01-02T10:29:24.327858Z","iopub.status.idle":"2025-01-02T10:31:56.311548Z","shell.execute_reply.started":"2025-01-02T10:29:24.327811Z","shell.execute_reply":"2025-01-02T10:31:56.310268Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Plot predicted scores against true scores with the base thresholds for converting score to sii\nsns.set_theme(style=\"white\")\nfig, axes = plt.subplots(1, 5, figsize=(14, 6))\n\nscatter1 = axes[0].scatter(train['PCIAT-PCIAT_Total'], oof_lgb, c=train[\"sii\"], cmap=\"autumn\", alpha=0.5)\naxes[0].set_xlabel(\"True Score\")\naxes[0].set_ylabel(\"OOF Predictions - LGBM\")\naxes[0].set_ylim(0,np.max(train['PCIAT-PCIAT_Total']))\naxes[0].set_xlim(0,np.max(train['PCIAT-PCIAT_Total']))\naxes[0].set_aspect('equal', adjustable='box')\n\nthresholds = [30, 50, 80]\nfor threshold in thresholds:\n    axes[0].axhline(threshold, color=\"blue\", linestyle=\"--\", lw=1)\n    axes[0].axvline(threshold, color=\"blue\", linestyle=\"--\", lw=1)\n\nscatter2 = axes[1].scatter(train['PCIAT-PCIAT_Total'], oof_xgb, c=train[\"sii\"], cmap=\"autumn\", alpha=0.5)\naxes[1].set_xlabel(\"True Score\")\naxes[1].set_ylabel(\"OOF Predictions - XGB\")\naxes[1].set_ylim(0,np.max(train['PCIAT-PCIAT_Total']))\naxes[1].set_xlim(0,np.max(train['PCIAT-PCIAT_Total']))\naxes[1].set_aspect('equal', adjustable='box')\n\nfor threshold in thresholds:\n    axes[1].axhline(threshold, color=\"blue\", linestyle=\"--\", lw=1)\n    axes[1].axvline(threshold, color=\"blue\", linestyle=\"--\", lw=1)\n\nscatter2 = axes[2].scatter(train['PCIAT-PCIAT_Total'], oof_xgb_2, c=train[\"sii\"], cmap=\"autumn\", alpha=0.5)\naxes[2].set_xlabel(\"True Score\")\naxes[2].set_ylabel(\"OOF Predictions - XGB (2)\")\naxes[2].set_ylim(0,np.max(train['PCIAT-PCIAT_Total']))\naxes[2].set_xlim(0,np.max(train['PCIAT-PCIAT_Total']))\naxes[2].set_aspect('equal', adjustable='box')\n\nfor threshold in thresholds:\n    axes[2].axhline(threshold, color=\"blue\", linestyle=\"--\", lw=1)\n    axes[2].axvline(threshold, color=\"blue\", linestyle=\"--\", lw=1)\n    \nscatter3 = axes[3].scatter(train['PCIAT-PCIAT_Total'], oof_cat, c=train[\"sii\"], cmap=\"autumn\", alpha=0.5)\naxes[3].set_xlabel(\"True Score\")\naxes[3].set_ylabel(\"OOF Predictions - Cat\")\naxes[3].set_ylim(0,np.max(train['PCIAT-PCIAT_Total']))\naxes[3].set_xlim(0,np.max(train['PCIAT-PCIAT_Total']))\naxes[3].set_aspect('equal', adjustable='box')\n\nfor threshold in thresholds:\n    axes[3].axhline(threshold, color=\"blue\", linestyle=\"--\", lw=1)\n    axes[3].axvline(threshold, color=\"blue\", linestyle=\"--\", lw=1)\n    \nscatter3 = axes[4].scatter(train['PCIAT-PCIAT_Total'], oof_xtrees, c=train[\"sii\"], cmap=\"autumn\", alpha=0.5)\naxes[4].set_xlabel(\"True Score\")\naxes[4].set_ylabel(\"OOF Predictions - ExtraTrees\")\naxes[4].set_ylim(0,np.max(train['PCIAT-PCIAT_Total']))\naxes[4].set_xlim(0,np.max(train['PCIAT-PCIAT_Total']))\naxes[4].set_aspect('equal', adjustable='box')\n\nfor threshold in thresholds:\n    axes[4].axhline(threshold, color=\"blue\", linestyle=\"--\", lw=1)\n    axes[4].axvline(threshold, color=\"blue\", linestyle=\"--\", lw=1)\n\nplt.tight_layout()\nplt.show()\n\nfig, axes = plt.subplots(1,3, figsize=(20,6))\n\nmodel_preds = pd.DataFrame({\n    'lgb': oof_lgb,\n    'xgb': oof_xgb,\n    'xgb2': oof_xgb_2,\n    'cat': oof_cat,\n    'xtrees': oof_xtrees\n})\ncorr_df = model_preds.corr()\nsns.heatmap(corr_df, annot=True, cmap=\"autumn\", cbar=False, linewidths=0.5, linecolor='black', ax=axes[0])\naxes[0].set_title(\"Correlation Between Models\")\n\naxes[1].scatter(\n    train['PCIAT-PCIAT_Total'], np.average(model_preds, axis=1, weights=[0.2, 0.2, 0.3, 0.1, 0.2]), c=train[\"sii\"], cmap=\"autumn\", alpha=0.5\n)\naxes[1].set_xlabel(\"True Score\")\naxes[1].set_ylabel(\"OOF Predictions - Ensemble\")\naxes[1].set_ylim(0,np.max(train['PCIAT-PCIAT_Total']))\naxes[1].set_xlim(0,np.max(train['PCIAT-PCIAT_Total']))\naxes[1].set_aspect('equal', adjustable='box')\nfor threshold in thresholds:\n    axes[1].axhline(threshold, color=\"blue\", linestyle=\"--\", lw=1)\n    axes[1].axvline(threshold, color=\"blue\", linestyle=\"--\", lw=1)\n\n\nmeta_model = LinearRegression(fit_intercept = False)\nmeta_model.fit(model_preds, train['PCIAT-PCIAT_Total'])\nmeta_model.coef_\nprint(f\"Coefs 'Meta Model': {meta_model.coef_}\")\n\nlgb_thresholds_ens = np.mean(np.array(lgb_thresholds), axis=0)\nxgb_thresholds_ens = np.mean(np.array(xgb_thresholds), axis=0)\nxgb_2_thresholds_ens = np.mean(np.array(xgb_2_thresholds), axis=0)\ncat_thresholds_ens = np.mean(np.array(cat_thresholds), axis=0)\nxtrees_thresholds_ens = np.mean(np.array(xtrees_thresholds), axis=0)\n\nthresholds_df = pd.DataFrame({\n    \"LGB Thresholds\": lgb_thresholds_ens,\n    \"XGB Thresholds\": xgb_thresholds_ens,\n    \"XGB 2 Thresholds\": xgb_2_thresholds_ens,\n    \"Cat Thresholds\": cat_thresholds_ens,\n    \"Xtrees Thresholds\": xtrees_thresholds_ens\n})\n\nsns.heatmap(thresholds_df, annot=True, cmap=\"autumn\", cbar=False, linewidths=0.5, linecolor='black', ax=axes[2])\naxes[2].set_title(\"Ensemble Thresholds Derived from CV\")\naxes[2].set_xticklabels(thresholds_df.columns, rotation=45)  \nplt.show()","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:31:56.312788Z","iopub.execute_input":"2025-01-02T10:31:56.313283Z","iopub.status.idle":"2025-01-02T10:31:58.919098Z","shell.execute_reply.started":"2025-01-02T10:31:56.313249Z","shell.execute_reply":"2025-01-02T10:31:58.917686Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# # Showing how sensitive QWK is with changes in threshold\n# scores = []\n# ts = []\n# m = 30\n# for i in tqdm(np.linspace(-1,10, 100)):\n#     thresholds = [m+i, 50, 80]\n#     pred = round_with_thresholds(oof_xgb, thresholds)\n#     score = cohen_kappa_score(train[\"sii\"], pred, weights='quadratic')\n#     ts.append(m+i)\n#     scores.append(score)\n    \n# plt.plot(ts, scores)\n# plt.title(\"Demonstration of QWK sensitivity to changes in threshold 0\")\n# plt.show()","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:31:58.920429Z","iopub.execute_input":"2025-01-02T10:31:58.920878Z","iopub.status.idle":"2025-01-02T10:31:59.574476Z","shell.execute_reply.started":"2025-01-02T10:31:58.920837Z","shell.execute_reply":"2025-01-02T10:31:59.573174Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"lgb_importances = pd.DataFrame({\n    'Feature': lgb_features,\n    'Importance': lgb_model.feature_importances_\n}).sort_values(by='Importance', ascending=False)\n\nxgb_importances = pd.DataFrame({\n    'Feature': xgb_features,\n    'Importance': xgb_model.feature_importances_\n}).sort_values(by='Importance', ascending=False)\n\ncat_importances = pd.DataFrame({\n    'Feature': cat_features,\n    'Importance': cat_model.feature_importances_\n}).sort_values(by='Importance', ascending=False)\n\n# Set the number of features to display\nn_top_features = 40\n\nfig, axes = plt.subplots(1, 3, figsize=(18, 8))\nsns.set_theme(style=\"whitegrid\")\n\nsns.barplot(ax=axes[0], data=lgb_importances.head(n_top_features),\n            x='Importance', y='Feature', palette=\"autumn\")\naxes[0].set_title('LightGBM Top Feature Importances')\n\nsns.barplot(ax=axes[1], data=xgb_importances.head(n_top_features),\n            x='Importance', y='Feature', palette=\"autumn\")\naxes[1].set_title('XGBoost Top Feature Importances')\n\nsns.barplot(ax=axes[2], data=cat_importances.head(n_top_features),\n            x='Importance', y='Feature', palette=\"autumn\")\naxes[2].set_title('CatBoost Top Feature Importances')\n\nplt.tight_layout()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:31:59.576082Z","iopub.execute_input":"2025-01-02T10:31:59.576508Z","iopub.status.idle":"2025-01-02T10:32:01.665485Z","shell.execute_reply.started":"2025-01-02T10:31:59.576477Z","shell.execute_reply":"2025-01-02T10:32:01.664296Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Features of low importance in XGB (1), LGBM and Cat\n(set(lgb_importances[-10:][\"Feature\"]) & set(xgb_importances[-10:][\"Feature\"]) & set(cat_importances[-10:][\"Feature\"]))","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:32:01.666614Z","iopub.execute_input":"2025-01-02T10:32:01.666952Z","iopub.status.idle":"2025-01-02T10:32:01.675312Z","shell.execute_reply.started":"2025-01-02T10:32:01.666915Z","shell.execute_reply":"2025-01-02T10:32:01.673583Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Apply the optimized thresholds to OOF predictions\noof_lgb = round_with_thresholds(oof_lgb, lgb_thresholds_ens)\nprint(f\"LGBM optimized Kappa: {cohen_kappa_score(train['sii'], oof_lgb, weights='quadratic')}\")\noof_xgb = round_with_thresholds(oof_xgb, xgb_thresholds_ens)\nprint(f\"XGB optimized Kappa: {cohen_kappa_score(train['sii'], oof_xgb, weights='quadratic')}\")\noof_xgb_2 = round_with_thresholds(oof_xgb_2, xgb_2_thresholds_ens)\nprint(f\"XGB2 optimized Kappa: {cohen_kappa_score(train['sii'], oof_xgb_2, weights='quadratic')}\")\noof_cat = round_with_thresholds(oof_cat, cat_thresholds_ens)\nprint(f\"CAT optimized Kappa: {cohen_kappa_score(train['sii'], oof_cat, weights='quadratic')}\")\noof_xtrees = round_with_thresholds(oof_xtrees, xtrees_thresholds_ens)\nprint(f\"ExtraTrees optimized Kappa: {cohen_kappa_score(train['sii'], oof_xtrees, weights='quadratic')}\")\n\nif voting:\n    oof_preds = np.array([oof_lgb, oof_xgb, oof_xgb_2, oof_cat, oof_xtrees])\n    voted_oof = stats.mode(oof_preds, axis=0).mode.flatten().astype(int)\n    final_oof = voted_oof\nelse: \n    weights = [0.2, 0.2, 0.3, 0.1, 0.2]\n    oof_preds = np.array([oof_lgb, oof_xgb, oof_xgb_2, oof_cat, oof_xtrees])\n    weighted_oof = np.average(oof_preds, axis=0, weights=weights)\n    final_oof = np.round(weighted_oof).astype(int)\n\n# Calculate Kappa score for voted OOF predictions\nkappa_score = cohen_kappa_score(train[\"sii\"], final_oof, weights='quadratic')\nprint(f\"Ensemble Kappa score: {kappa_score}\")\n# Plot confusion matrix\nconf_matrix = confusion_matrix(train[\"sii\"], final_oof)\nsns.set_theme(style=\"whitegrid\")\nplt.figure(figsize=(8, 6))\nsns.heatmap(conf_matrix, annot=True, fmt=\"d\", cmap=\"autumn\", cbar=False, linewidths=0.5, linecolor='black')\nplt.title('Confusion Matrix', fontsize=16)\nplt.xlabel('Predicted', fontsize=12)\nplt.ylabel('True', fontsize=12)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:32:01.676954Z","iopub.execute_input":"2025-01-02T10:32:01.677307Z","iopub.status.idle":"2025-01-02T10:32:02.150255Z","shell.execute_reply.started":"2025-01-02T10:32:01.677275Z","shell.execute_reply":"2025-01-02T10:32:02.149178Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Apply the optimized thresholds to test predictions\ntest_lgb = round_with_thresholds(test_lgb, lgb_thresholds_ens)\ntest_xgb = round_with_thresholds(test_xgb, xgb_thresholds_ens)\ntest_xgb_2 = round_with_thresholds(test_xgb_2, xgb_2_thresholds_ens)\ntest_cat = round_with_thresholds(test_cat, cat_thresholds_ens)\ntest_xtrees = round_with_thresholds(test_xtrees, xtrees_thresholds_ens)\n\nif voting:\n    test_preds = np.array([test_lgb, test_xgb, test_xgb_2, test_cat, test_xtrees])\n    voted_test = stats.mode(test_preds, axis=0).mode.flatten().astype(int)\n    final_test = np.round(voted_test).astype(int)\nelse:\n    test_preds = np.array([test_lgb, test_xgb, test_xgb_2, test_cat, test_xtrees])\n    weighted_test = np.average(test_preds, axis=0, weights=weights)\n    final_test = np.round(weighted_test).astype(int)\n\nsubmission = pd.read_csv(SUBMISSION_PATH)\nsubmission['sii'] = final_test\n\n# Define the output path for the submission file\nsubmission_path = os.path.join(OUTPUT_PATH, \"submission.csv\")\n\n# Save the submission file\nsubmission.to_csv(submission_path, index=False)","metadata":{"execution":{"iopub.status.busy":"2025-01-02T10:32:02.151369Z","iopub.execute_input":"2025-01-02T10:32:02.151759Z","iopub.status.idle":"2025-01-02T10:32:02.177778Z","shell.execute_reply.started":"2025-01-02T10:32:02.151720Z","shell.execute_reply":"2025-01-02T10:32:02.176437Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}