{"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":"nvidiaTeslaT4","dataSources":[{"sourceId":81933,"databundleVersionId":9643020,"sourceType":"competition"},{"sourceId":197506774,"sourceType":"kernelVersion"}],"dockerImageVersionId":30823,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"### This notebook is fork of [this](https://www.kaggle.com/code/tubotubo/starter-notebook-multi-target-prediction).\n\nThe difference is stats features from parquets. I used **polars** with **GPU** for feature engineering.\n\nThis notebook illustrates how to use polars-gpu briefly.","metadata":{}},{"cell_type":"code","source":"!nvidia-smi","metadata":{"execution":{"iopub.status.busy":"2024-12-20T08:14:15.390367Z","iopub.execute_input":"2024-12-20T08:14:15.390724Z","iopub.status.idle":"2024-12-20T08:14:15.620229Z","shell.execute_reply.started":"2024-12-20T08:14:15.390692Z","shell.execute_reply":"2024-12-20T08:14:15.619329Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"!pip install /kaggle/input/polars-gpu-1-7-1/cupy_cuda12x-13.3.0-cp310-cp310-manylinux2014_x86_64.whl\n!pip install /kaggle/input/polars-gpu-1-7-1/rmm_cu12-24.8.2-cp310-cp310-manylinux_2_24_x86_64.manylinux_2_28_x86_64.whl\n!pip install /kaggle/input/polars-gpu-1-7-1/cudf_cu12-24.8.3-cp310-cp310-manylinux_2_28_x86_64.whl\n!pip install /kaggle/input/polars-gpu-1-7-1/polars-1.7.1-cp38-abi3-manylinux_2_17_x86_64.manylinux2014_x86_64.whl\n!pip install /kaggle/input/polars-gpu-1-7-1/cudf_polars_cu12-24.8.3-py3-none-any.whl\n!pip install cudf-polars","metadata":{"execution":{"iopub.status.busy":"2024-12-20T08:14:15.621293Z","iopub.execute_input":"2024-12-20T08:14:15.621547Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"!pip install /kaggle/input/polars-gpu-1-7-1/*.whl","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import warnings\nfrom functools import partial\nfrom pathlib import Path\n\nfrom tqdm import tqdm\nimport matplotlib.pyplot as plt\nimport numpy as np\nimport optuna\nimport polars as pl\nimport polars.selectors as cs\nfrom catboost import CatBoostRegressor, MultiTargetCustomMetric\nfrom numpy.typing import ArrayLike, NDArray\nfrom polars.testing import assert_frame_equal\nfrom sklearn.base import BaseEstimator\nfrom sklearn.metrics import cohen_kappa_score\nfrom sklearn.model_selection import StratifiedKFold\n\nwarnings.filterwarnings(\"ignore\", message=\"Failed to optimize method\")\n\nDATA_DIR = Path(\"/kaggle/input/child-mind-institute-problematic-internet-use\")\nTARGET_COLS = [\n    \"PCIAT-PCIAT_01\",\n    \"PCIAT-PCIAT_02\",\n    \"PCIAT-PCIAT_03\",\n    \"PCIAT-PCIAT_04\",\n    \"PCIAT-PCIAT_05\",\n    \"PCIAT-PCIAT_06\",\n    \"PCIAT-PCIAT_07\",\n    \"PCIAT-PCIAT_08\",\n    \"PCIAT-PCIAT_09\",\n    \"PCIAT-PCIAT_10\",\n    \"PCIAT-PCIAT_11\",\n    \"PCIAT-PCIAT_12\",\n    \"PCIAT-PCIAT_13\",\n    \"PCIAT-PCIAT_14\",\n    \"PCIAT-PCIAT_15\",\n    \"PCIAT-PCIAT_16\",\n    \"PCIAT-PCIAT_17\",\n    \"PCIAT-PCIAT_18\",\n    \"PCIAT-PCIAT_19\",\n    \"PCIAT-PCIAT_20\",\n    \"PCIAT-PCIAT_Total\",\n    \"sii\",\n]\n\nFEATURE_COLS = [\n    \"Basic_Demos-Enroll_Season\",\n    \"Basic_Demos-Age\",\n    \"Basic_Demos-Sex\",\n    \"CGAS-Season\",\n    \"CGAS-CGAS_Score\",\n    \"Physical-Season\",\n    \"Physical-BMI\",\n    \"Physical-Height\",\n    \"Physical-Weight\",\n    \"Physical-Waist_Circumference\",\n    \"Physical-Diastolic_BP\",\n    \"Physical-HeartRate\",\n    \"Physical-Systolic_BP\",\n    \"Fitness_Endurance-Season\",\n    \"Fitness_Endurance-Max_Stage\",\n    \"Fitness_Endurance-Time_Mins\",\n    \"Fitness_Endurance-Time_Sec\",\n    \"FGC-Season\",\n    \"FGC-FGC_CU\",\n    \"FGC-FGC_CU_Zone\",\n    \"FGC-FGC_GSND\",\n    \"FGC-FGC_GSND_Zone\",\n    \"FGC-FGC_GSD\",\n    \"FGC-FGC_GSD_Zone\",\n    \"FGC-FGC_PU\",\n    \"FGC-FGC_PU_Zone\",\n    \"FGC-FGC_SRL\",\n    \"FGC-FGC_SRL_Zone\",\n    \"FGC-FGC_SRR\",\n    \"FGC-FGC_SRR_Zone\",\n    \"FGC-FGC_TL\",\n    \"FGC-FGC_TL_Zone\",\n    \"BIA-Season\",\n    \"BIA-BIA_Activity_Level_num\",\n    \"BIA-BIA_BMC\",\n    \"BIA-BIA_BMI\",\n    \"BIA-BIA_BMR\",\n    \"BIA-BIA_DEE\",\n    \"BIA-BIA_ECW\",\n    \"BIA-BIA_FFM\",\n    \"BIA-BIA_FFMI\",\n    \"BIA-BIA_FMI\",\n    \"BIA-BIA_Fat\",\n    \"BIA-BIA_Frame_num\",\n    \"BIA-BIA_ICW\",\n    \"BIA-BIA_LDM\",\n    \"BIA-BIA_LST\",\n    \"BIA-BIA_SMM\",\n    \"BIA-BIA_TBW\",\n    \"PAQ_A-Season\",\n    \"PAQ_A-PAQ_A_Total\",\n    \"PAQ_C-Season\",\n    \"PAQ_C-PAQ_C_Total\",\n    \"SDS-Season\",\n    \"SDS-SDS_Total_Raw\",\n    \"SDS-SDS_Total_T\",\n    \"PreInt_EduHx-Season\",\n    \"PreInt_EduHx-computerinternet_hoursday\",\n    \n    # stats features from parquets\n    \"X_min\",\n    \"Y_min\",\n    \"Z_min\",\n    \"enmo_min\",\n    \"anglez_min\",\n    \"light_min\",\n    \"battery_voltage_min\",\n    \"X_mean\",\n    \"Y_mean\",\n    \"Z_mean\",\n    \"enmo_mean\",\n    \"anglez_mean\",\n    \"light_mean\",\n    \"battery_voltage_mean\",\n    \"X_max\",\n    \"Y_max\",\n    \"Z_max\",\n    \"enmo_max\",\n    \"anglez_max\",\n    \"light_max\",\n    \"battery_voltage_max\",\n    \"X_std\",\n    \"Y_std\",\n    \"Z_std\",\n    \"enmo_std\",\n    \"anglez_std\",\n    \"light_std\",\n    \"battery_voltage_std\",\n]","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Load data\ntrain = pl.read_csv(DATA_DIR / \"train.csv\")\ntest = pl.read_csv(DATA_DIR / \"test.csv\")\ntrain_test = pl.concat([train, test], how=\"diagonal\")\n\nIS_TEST = test.height <= 100\n\nassert_frame_equal(train, train_test[: train.height].select(train.columns))\nassert_frame_equal(test, train_test[train.height :].select(test.columns))","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Cast string columns to categorical\ntrain_test = train_test.with_columns(cs.string().cast(pl.Categorical).fill_null(\"NAN\"))\ntrain = train_test[: train.height]\ntest = train_test[train.height :]","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def split_array(ar, n_group):\n    for i_chunk in range(n_group):\n        yield ar[i_chunk * len(ar) // n_group : (i_chunk + 1) * len(ar) // n_group]\n\ndef agg_parquets(files):\n    cols = [\"X\", \"Y\", \"Z\", \"enmo\", \"anglez\", \"light\", \"battery_voltage\"]\n    aggs = []\n    files_chunks = list(split_array(files, 10))\n    for files_tmp in tqdm(files_chunks):\n        if len(files_tmp) == 0:\n            continue\n        dfs = []\n        for file in files_tmp:\n            df = pl.scan_parquet(file)\n            df = df.with_columns(pl.lit(file.parts[-1].split(\"=\")[1]).alias(\"id\"))\n            dfs.append(df)\n        df = pl.concat(dfs)\n        agg = (\n            df.group_by(\"id\")\n            .agg(\n                [pl.col(c).cast(pl.Float32).min().alias(f\"{c}_min\") for c in cols]\n                + [pl.col(c).cast(pl.Float32).mean().alias(f\"{c}_mean\") for c in cols]\n                + [pl.col(c).cast(pl.Float32).max().alias(f\"{c}_max\") for c in cols]\n                + [pl.col(c).cast(pl.Float32).std().alias(f\"{c}_std\") for c in cols]\n            )\n            .collect(engine=\"gpu\")\n        )\n        aggs.append(agg)\n    return pl.concat(aggs)\n\ntrain_agg = agg_parquets(sorted(Path(\"/kaggle/input/child-mind-institute-problematic-internet-use/series_train.parquet\").glob(\"*\")))\ntest_agg = agg_parquets(sorted(Path(\"/kaggle/input/child-mind-institute-problematic-internet-use/series_test.parquet\").glob(\"*\")))","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"train_agg","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"test_agg","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"train = train.join(train_agg.with_columns(pl.col(\"id\").cast(pl.Categorical)), on='id', how='left')\ntest = test.join(test_agg.with_columns(pl.col(\"id\").cast(pl.Categorical)), on='id', how='left')","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"test","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ignore rows with null values in TARGET_COLS\ntrain_without_null = train.drop_nulls(subset=TARGET_COLS)\nX = train_without_null.select(FEATURE_COLS)\nX_test = test.select(FEATURE_COLS)\ny = train_without_null.select(TARGET_COLS)\ny_sii = y.get_column(\"sii\").to_numpy()  # ground truth\ncat_features = X.select(cs.categorical()).columns","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"X","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"X_test","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"class MultiTargetQWK(MultiTargetCustomMetric):\n    def get_final_error(self, error, weight):\n        return np.sum(error)  # / np.sum(weight)\n\n    def is_max_optimal(self):\n        # if True, the bigger the better\n        return True\n\n    def evaluate(self, approxes, targets, weight):\n        # approxes: 予測値 (shape: [ターゲット数, サンプル数])\n        # targets: 実際の値 (shape: [ターゲット数, サンプル数])\n        # weight: サンプルごとの重み (Noneも可)\n\n        approx = np.clip(approxes[-1], 0, 3).round().astype(int)\n        target = targets[-1]\n\n        qwk = cohen_kappa_score(target, approx, weights=\"quadratic\")\n\n        return qwk, 1\n\n    def get_custom_metric_name(self):\n        return \"MultiTargetQWK\"\n\n\nclass OptimizedRounder:\n    \"\"\"\n    A class for optimizing the rounding of continuous predictions into discrete class labels using Optuna.\n    The optimization process maximizes the Quadratic Weighted Kappa score by learning thresholds that separate\n    continuous predictions into class intervals.\n\n    Args:\n        n_classes (int): The number of discrete class labels.\n        n_trials (int, optional): The number of trials for the Optuna optimization. Defaults to 100.\n\n    Attributes:\n        n_classes (int): The number of discrete class labels.\n        labels (NDArray[np.int_]): An array of class labels from 0 to `n_classes - 1`.\n        n_trials (int): The number of optimization trials.\n        metric (Callable): The Quadratic Weighted Kappa score metric used for optimization.\n        thresholds (List[float]): The optimized thresholds learned after calling `fit()`.\n\n    Methods:\n        fit(y_pred: NDArray[np.float_], y_true: NDArray[np.int_]) -> None:\n            Fits the rounding thresholds based on continuous predictions and ground truth labels.\n\n            Args:\n                y_pred (NDArray[np.float_]): Continuous predictions that need to be rounded.\n                y_true (NDArray[np.int_]): Ground truth class labels.\n\n            Returns:\n                None\n\n        predict(y_pred: NDArray[np.float_]) -> NDArray[np.int_]:\n            Predicts discrete class labels by rounding continuous predictions using the fitted thresholds.\n            `fit()` must be called before `predict()`.\n\n            Args:\n                y_pred (NDArray[np.float_]): Continuous predictions to be rounded.\n\n            Returns:\n                NDArray[np.int_]: Predicted class labels.\n\n        _normalize(y: NDArray[np.float_]) -> NDArray[np.float_]:\n            Normalizes the continuous values to the range [0, `n_classes - 1`].\n\n            Args:\n                y (NDArray[np.float_]): Continuous values to be normalized.\n\n            Returns:\n                NDArray[np.float_]: Normalized values.\n\n    References:\n        - This implementation uses Optuna for threshold optimization.\n        - Quadratic Weighted Kappa is used as the evaluation metric.\n    \"\"\"\n\n    def __init__(self, n_classes: int, n_trials: int = 100):\n        self.n_classes = n_classes\n        self.labels = np.arange(n_classes)\n        self.n_trials = n_trials\n        self.metric = partial(cohen_kappa_score, weights=\"quadratic\")\n\n    def fit(self, y_pred: NDArray[np.float_], y_true: NDArray[np.int_]) -> None:\n        y_pred = self._normalize(y_pred)\n\n        def objective(trial: optuna.Trial) -> float:\n            thresholds = []\n            for i in range(self.n_classes - 1):\n                low = max(thresholds) if i > 0 else min(self.labels)\n                high = max(self.labels)\n                th = trial.suggest_float(f\"threshold_{i}\", low, high)\n                thresholds.append(th)\n            try:\n                y_pred_rounded = np.digitize(y_pred, thresholds)\n            except ValueError:\n                return -100\n            return self.metric(y_true, y_pred_rounded)\n\n        optuna.logging.disable_default_handler()\n        study = optuna.create_study(direction=\"maximize\")\n        study.optimize(\n            objective,\n            n_trials=self.n_trials,\n        )\n        self.thresholds = [study.best_params[f\"threshold_{i}\"] for i in range(self.n_classes - 1)]\n\n    def predict(self, y_pred: NDArray[np.float_]) -> NDArray[np.int_]:\n        assert hasattr(self, \"thresholds\"), \"fit() must be called before predict()\"\n        y_pred = self._normalize(y_pred)\n        return np.digitize(y_pred, self.thresholds)\n\n    def _normalize(self, y: NDArray[np.float_]) -> NDArray[np.float_]:\n        # normalize y_pred to [0, n_classes - 1]\n        return (y - y.min()) / (y.max() - y.min()) * (self.n_classes - 1)","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# setting catboost parameters\nparams = dict(\n    loss_function=\"MultiRMSE\",\n    eval_metric=MultiTargetQWK(),\n    iterations=1 if IS_TEST else 100000,\n    learning_rate=0.1,\n    depth=5,\n    early_stopping_rounds=50,\n)\n\n# Cross-validation\nskf = StratifiedKFold(n_splits=95, shuffle=True, random_state=52)\nmodels: list[CatBoostRegressor] = []\ny_pred = np.full((X.height, len(TARGET_COLS)), fill_value=np.nan)\nfor train_idx, val_idx in skf.split(X, y_sii):\n    X_train: pl.DataFrame\n    X_val: pl.DataFrame\n    y_train: pl.DataFrame\n    y_val: pl.DataFrame\n    X_train, X_val = X[train_idx], X[val_idx]\n    y_train, y_val = y[train_idx], y[val_idx]\n\n    # train model\n    model = CatBoostRegressor(**params)\n    model.fit(\n        X_train.to_pandas(),\n        y_train.to_pandas(),\n        eval_set=(X_val.to_pandas(), y_val.to_pandas()),\n        cat_features=cat_features,\n        verbose=False,\n    )\n    models.append(model)\n\n    # predict\n    y_pred[val_idx] = model.predict(X_val.to_pandas())\n\nassert np.isnan(y_pred).sum() == 0\n# Optimize thresholds\noptimizer = OptimizedRounder(n_classes=4, n_trials=300)\ny_pred_total = y_pred[:, TARGET_COLS.index(\"PCIAT-PCIAT_Total\")]\noptimizer.fit(y_pred_total, y_sii)\ny_pred_rounded = optimizer.predict(y_pred_total)\n\n# Calculate QWK\nqwk = cohen_kappa_score(y_sii, y_pred_rounded, weights=\"quadratic\")\nprint(f\"Cross-Validated QWK Score: {qwk}\")","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"feature_importance = np.mean([model.get_feature_importance() for model in models], axis=0)\nsorted_idx = np.argsort(feature_importance)\nsorted_idx = sorted_idx[-30:]\nfig = plt.figure(figsize=(12, 10))\nplt.barh(range(len(sorted_idx)), feature_importance[sorted_idx], align=\"center\")\nplt.yticks(range(len(sorted_idx)), np.array(X_test.columns)[sorted_idx])\nplt.title(\"Feature Importance\")","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"class AvgModel:\n    def __init__(self, models: list[BaseEstimator]):\n        self.models = models\n\n    def predict(self, X: ArrayLike) -> NDArray[np.int_]:\n        preds: list[NDArray[np.int_]] = []\n        for model in self.models:\n            pred = model.predict(X)\n            preds.append(pred)\n\n        return np.mean(preds, axis=0)","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"avg_model = AvgModel(models)\ntest_pred = avg_model.predict(X_test.to_pandas())[:, TARGET_COLS.index(\"PCIAT-PCIAT_Total\")]\ntest_pred_rounded = optimizer.predict(test_pred)\ntest.select(\"id\").with_columns(\n    pl.Series(\"sii\", pl.Series(\"sii\", test_pred_rounded)),\n).write_csv(\"submission.csv\")","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}