{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import os\n\nimport numpy as np\nimport pandas as pd\nimport seaborn as sns\nimport matplotlib.pyplot as plt\n\nfrom typing import List, Tuple, Optional\nfrom tqdm import tqdm\nfrom datetime import datetime\n\nfrom IPython.core.interactiveshell import InteractiveShell\n\nfrom pathlib import Path\n\nimport polars as pl\n\nInteractiveShell.ast_node_interactivity = \"all\"","metadata":{"execution":{"iopub.status.busy":"2023-11-13T16:49:34.907487Z","iopub.execute_input":"2023-11-13T16:49:34.907982Z","iopub.status.idle":"2023-11-13T16:49:34.916573Z","shell.execute_reply.started":"2023-11-13T16:49:34.907938Z","shell.execute_reply":"2023-11-13T16:49:34.915123Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"1st stage (Predict whether any event occurs within a window)\n- The training input should only select the window near the step the events occur. Limit the snapshots step_event +/- training_window as downsampling. Self-joining on Polars is too expensive.\n- Prediction: [0.5 0.3 0.2 0.1 0.75] ; each point refers to the window of size = t\n- Classifying window strategy: Either using the threshold or the biggest gap in prediction\n- Denote 1 = awake, 0 = onset; Diff(Pred) > 0 (awake) otherwise (onset)\n\n2nd stage (Predict when the event occurs within a window)\n- Input data format:\nstep_in_window   Event  Feature (Reused)\n    1               0      ...\n    2               0      ...\n    3               1      ...\n    4               0      ...\n- Target distribution: discrete\n- Example features: normalized features within a series or window, any feature value relative to starting point of a window (Recency)\n- 2 approaches: Binary Classification, Ordinal Regression\n\nResult\n- Predict \"onset\" and \"wakeup\" at t and t+63 with the probability = P(event occurs | event in window) * P(event in window)\n- P(event in window) is normalized by softmax","metadata":{}},{"cell_type":"code","source":"dt_transforms = [\n    pl.col('timestamp').str.to_datetime(), \n    (pl.col('timestamp').str.to_datetime().dt.year()-2000).cast(pl.UInt8).alias('year'), \n    pl.col('timestamp').str.to_datetime().dt.month().cast(pl.UInt8).alias('month'),\n    pl.col('timestamp').str.to_datetime().dt.day().cast(pl.UInt8).alias('day'), \n    pl.col('timestamp').str.to_datetime().dt.hour().cast(pl.UInt8).alias('hour'),\n    pl.col('timestamp').str.to_datetime().dt.minute().cast(pl.UInt8).alias('minute')\n]\n\ndata_transforms = [\n    pl.col('anglez').cast(pl.Int16), # Casting anglez to 16 bit integer\n    (pl.col('enmo')*1000).cast(pl.UInt16), # Convert enmo to 16 bit uint\n]\n\ntrain_series_pl = pl.scan_parquet('/kaggle/input/child-mind-institute-detect-sleep-states/train_series.parquet').with_columns(dt_transforms + data_transforms)\ntrain_events_pl = pl.read_csv('/kaggle/input/child-mind-institute-detect-sleep-states/train_events.csv').with_columns(dt_transforms).filter(pl.col(\"step\").is_null())","metadata":{"execution":{"iopub.status.busy":"2023-11-13T16:49:36.847737Z","iopub.execute_input":"2023-11-13T16:49:36.848182Z","iopub.status.idle":"2023-11-13T16:49:36.916501Z","shell.execute_reply.started":"2023-11-13T16:49:36.848134Z","shell.execute_reply":"2023-11-13T16:49:36.915384Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"FEATURE_NAMES = [\n    \"anglez\",\n    \"enmo\",\n    \"step\",\n    \"hour_sin\",\n    \"hour_cos\",\n    \"month_sin\",\n    \"month_cos\",\n    \"minute_sin\",\n    \"minute_cos\",\n    \"anglez_sin\",\n    \"anglez_cos\",\n]\n\ndef add_rolling_features(x: pl.Expr, rolling_steps: int, name: str) -> list[pl.Expr]:\n    # Center = True setting the labels at the center of the window, otherwise the window is forward-moving\n    \n    features = [x.rolling_mean(rolling_steps, center=True, min_periods=1).abs().cast(pl.Float32).alias(f'{name}_rolling_abs_mean')]\n    features += [x.rolling_max(rolling_steps, center=True, min_periods=1).abs().cast(pl.Float32).alias(f'{name}_rolling_abs_max')]\n    features += [x.rolling_min(rolling_steps, center=True, min_periods=1).abs().cast(pl.Float32).alias(f'{name}_rolling_abs_min')]\n    features += [x.diff().abs().rolling_std(rolling_steps, center=True, min_periods=1).cast(pl.Float32).alias(f'{name}_rolling_abs_diff_std')]\n    features += [x.diff().abs().rolling_mean(rolling_steps, center=True, min_periods=1).cast(pl.Float32).alias(f'{name}_rolling_abs_diff_mean')]\n\n    return features\n\n# def add_lag_features(x: pl.Expr)\n\n# Create UDF using pl.Expr: convert anglez to rad\ndef to_coord(x: pl.Expr, max_: int, name: str) -> list[pl.Expr]:\n    rad = 2 * np.pi * (x % max_) / max_\n    x_sin = rad.sin()\n    x_cos = rad.cos()\n    return [x_sin.alias(f\"{name}_sin\"), x_cos.alias(f\"{name}_cos\")]\n\ndef deg_to_rad(x: pl.Expr) -> pl.Expr:\n    return np.pi / 180 * x\n\n# def to_signal(hour: pl.Expr, minute: pl.Expr) -> list[pl.Expr]:\n#     signal_awake = dict(zip(range(1440), np.sin(np.linspace(0, np.pi, 1440) + 0.208 * np.pi) ** 24))\n#     signal_onset = dict(zip(range(1440), np.sin(np.linspace(0, np.pi, 1440) + 0.555 * np.pi) ** 24))\n#     return [(hour * 60 + minute).map(signal_awake).cast(pl.Float64).alias(\"signal_awake\"), (hour * 60 + minute).map(signal_onset).cast(pl.Float64).alias(\"signal_onset\")]\n\ndef add_feature(series_df: pl.DataFrame) -> pl.DataFrame:\n\n    signal_awake_pl = pl.from_records([list(range(1440)), list(np.sin(np.linspace(0, np.pi, 1440) + 0.208 * np.pi) ** 24)], schema = [\"hour_minute\", \"signal_awake\"])\n    signal_onset_pl = pl.from_records([list(range(1440)), list(np.sin(np.linspace(0, np.pi, 1440) + 0.555 * np.pi) ** 24)], schema = [\"hour_minute\", \"signal_onset\"])\n\n    series_df = (\n        series_df\n        # .with_row_count(\"step\")\n        .with_columns(deg_to_rad(pl.col(\"anglez\")).alias(\"anglez_rad\"))\n        .with_columns(\n            *to_coord(pl.col(\"timestamp\").dt.hour(), 24, \"hour\"),\n            *to_coord(pl.col(\"timestamp\").dt.month(), 12, \"month\"),\n            *to_coord(pl.col(\"timestamp\").dt.minute(), 60, \"minute\"),\n            # (pl.col(\"step\") / pl.count(\"step\")).alias(\"step_pct\"),\n            pl.col('anglez_rad').sin().alias('anglez_sin'),\n            pl.col('anglez_rad').cos().alias('anglez_cos')\n        )\n        .with_columns(((pl.col(\"enmo\") - pl.col(\"enmo\").mean().over([\"series_id\"])) / pl.col(\"enmo\").std().over([\"series_id\"])).alias(\"standardized_enmo\"))\n        .with_columns(\n            *add_rolling_features(pl.col(\"anglez\"), 60, \"anglez\"),\n            *add_rolling_features(pl.col(\"enmo\"), 60, \"enmo\"),\n            (pl.col(\"hour\")*60 + pl.col(\"minute\")).cast(pl.Int64).alias(\"hour_minute\"),\n            (pl.col(\"enmo\").diff() == 0).cast(pl.Int8).alias(\"is_enmo_diff_zero\"),\n            (pl.col(\"anglez\").diff() == 0).cast(pl.Int8).alias(\"is_anglez_diff_zero\")\n        )\n        .lazy().join(signal_awake_pl.select([\"hour_minute\", \"signal_awake\"]).lazy(), on = ['hour_minute'], how = \"left\")\n        .lazy().join(signal_onset_pl.select([\"hour_minute\", \"signal_onset\"]).lazy(), on = ['hour_minute'], how = \"left\")\n        .with_columns(\n            pl.col(\"is_enmo_diff_zero\").rolling_sum(120, center=True, min_periods=1).cast(pl.Float32).alias('rolling_sum_enmo_diff_zero'),\n            pl.col(\"is_anglez_diff_zero\").rolling_sum(120, center=True, min_periods=1).cast(pl.Float32).alias('rolling_sum_anglez_diff_zero')\n        )\n    )\n    return series_df","metadata":{"execution":{"iopub.status.busy":"2023-11-13T16:50:48.676512Z","iopub.execute_input":"2023-11-13T16:50:48.676954Z","iopub.status.idle":"2023-11-13T16:50:48.704796Z","shell.execute_reply.started":"2023-11-13T16:50:48.676918Z","shell.execute_reply":"2023-11-13T16:50:48.703115Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def sample_data_pl(series_pl: pl.DataFrame, event_series_pl: pl.DataFrame, n_samples: int):\n\n    sample_series_ids = (\n        series_pl\n        .select(\"series_id\")\n        .unique()\n        .collect(streaming=True)\n        .sample(n_samples)\n        .get_column(\"series_id\")\n    )\n    \n    return series_pl.filter(pl.col(\"series_id\").is_in(sample_series_ids)), event_series_pl.filter(pl.col(\"series_id\").is_in(sample_series_ids))\n\nsample_series_pl, sample_events_pl = sample_data_pl(train_series_pl, train_events_pl, 10)","metadata":{"execution":{"iopub.status.busy":"2023-11-13T16:50:25.488298Z","iopub.execute_input":"2023-11-13T16:50:25.488768Z","iopub.status.idle":"2023-11-13T16:50:31.067689Z","shell.execute_reply.started":"2023-11-13T16:50:25.488728Z","shell.execute_reply":"2023-11-13T16:50:31.066685Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sample_series_pl, sample_events_pl = sample_data_pl(train_series_pl, train_events_pl, 2)\nsample_series_pl_with_features = add_feature(sample_series_pl)\nsample_series_pl_with_features.collect()","metadata":{"execution":{"iopub.status.busy":"2023-11-13T16:50:51.336468Z","iopub.execute_input":"2023-11-13T16:50:51.337437Z","iopub.status.idle":"2023-11-13T16:51:09.756910Z","shell.execute_reply.started":"2023-11-13T16:50:51.337396Z","shell.execute_reply":"2023-11-13T16:51:09.755469Z"},"trusted":true},"execution_count":null,"outputs":[]}]}