{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"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":53666,"databundleVersionId":6589269,"sourceType":"competition"},{"sourceId":152619799,"sourceType":"kernelVersion"},{"sourceId":153811367,"sourceType":"kernelVersion"},{"sourceId":157569355,"sourceType":"kernelVersion"},{"sourceId":157604909,"sourceType":"kernelVersion"},{"sourceId":157869245,"sourceType":"kernelVersion"}],"dockerImageVersionId":30559,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# ライブラリ / 汎用関数","metadata":{}},{"cell_type":"code","source":"! pip install japanize_matplotlib --quiet","metadata":{"execution":{"iopub.status.busy":"2024-01-14T07:58:22.300457Z","iopub.execute_input":"2024-01-14T07:58:22.301527Z","iopub.status.idle":"2024-01-14T07:58:40.249900Z","shell.execute_reply.started":"2024-01-14T07:58:22.301471Z","shell.execute_reply":"2024-01-14T07:58:40.248496Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import pandas as pd\nimport pandas.tseries.offsets as offsets\nimport numpy as np\nimport pyarrow as pa\nimport pyarrow.parquet as pq\nimport matplotlib.pyplot as plt\nimport japanize_matplotlib\nimport seaborn as sns\nimport os\nimport gc\nimport time\nimport IPython.display\nimport collections\n\nimport tensorflow as tf\nfrom tensorflow.keras.models import Sequential\nfrom tensorflow.keras.layers import LSTM, Dense, Conv1D, Conv2D, Conv3D, Dropout\nfrom tensorflow.keras.utils import to_categorical, timeseries_dataset_from_array\nfrom tensorflow.keras.callbacks import ModelCheckpoint, EarlyStopping\nfrom tensorflow.keras.optimizers import Adam\nfrom tensorflow.keras.losses import BinaryCrossentropy\n\nfrom sklearn.metrics import confusion_matrix, accuracy_score, precision_score, recall_score, f1_score\nfrom sklearn.preprocessing import PolynomialFeatures\nfrom sklearn.linear_model import LinearRegression\n\nfrom scipy.optimize import curve_fit\nfrom scipy.stats import norm, ttest_ind\n\nfrom datetime import timedelta\nfrom scipy.stats import norm\nfrom pathlib import Path\nfrom tqdm import tqdm","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2024-01-14T07:58:40.252351Z","iopub.execute_input":"2024-01-14T07:58:40.252745Z","iopub.status.idle":"2024-01-14T07:58:57.761529Z","shell.execute_reply.started":"2024-01-14T07:58:40.252711Z","shell.execute_reply":"2024-01-14T07:58:57.760625Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def graphs(df, id_map, start=0, end=0):\n    if not end:\n        df = df[df['id_map']==id_map].copy().reset_index(drop=True)\n    else:\n        df = df[df['id_map']==id_map].copy().reset_index(drop=True).iloc[start:end]\n    steps = df[df['time'].apply(lambda x: str(x)[-8:]=='00:00:00')]['step']\n    flag = len(df) != df['step'].iat[-1]\n\n    _, axes = plt.subplots(2, 1, figsize=(20,6))\n    axes[0].set_title(id_map)\n    a = 'anglez' if 'anglez' in df.columns else 'a_b_var'\n    e = 'enmo' if 'enmo' in df.columns else 'e_b_mean'        \n    for ax, col in zip(axes.ravel(), [a, e]):\n        if 'state' in df.columns:\n            if col == 'anglez':\n                df_state = df['state'].apply(lambda x: 50 if x==1 else -50 if x==0 else 0)\n            else:\n                df_state = df['state'].apply(lambda x: 2 if x==1 else 0 if x==0 else 1)\n        X = df.index if flag else df['step']\n        ax.plot(X, df[col], linewidth=0.5, marker='', color='0.5', label=col)\n        if 'state' in df.columns:\n            ax.plot(X, df_state, linewidth=1, marker='', color='b', label='state')\n        for step in steps:\n            X = df[df['step']==step].index[0] if flag else step\n            ax.axvline(x=X, linewidth=1, color='r', alpha=0.5)\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2024-01-14T07:58:57.762839Z","iopub.execute_input":"2024-01-14T07:58:57.763801Z","iopub.status.idle":"2024-01-14T07:58:57.777535Z","shell.execute_reply.started":"2024-01-14T07:58:57.763753Z","shell.execute_reply":"2024-01-14T07:58:57.776350Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def time_to_sec(x):\n    return x.timestamp() % SEC_OF_DAY\n\ndef time_to_step(x):\n    return (x.timestamp() % SEC_OF_DAY) / 5","metadata":{"execution":{"iopub.status.busy":"2024-01-14T07:58:57.779910Z","iopub.execute_input":"2024-01-14T07:58:57.780248Z","iopub.status.idle":"2024-01-14T07:58:57.800044Z","shell.execute_reply.started":"2024-01-14T07:58:57.780221Z","shell.execute_reply":"2024-01-14T07:58:57.798874Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Const","metadata":{}},{"cell_type":"code","source":"HOME = Path('/kaggle')\nORIGIN = HOME / 'input/child-mind-institute-detect-sleep-states'\nCUSTOM = HOME / 'input/zzz-eda-dataset'\nWORKING = HOME / 'working'\n\nMODEL = []\nfor i in range(1,4):\n    path = f'input/zzz-model-{i}-ipynb'\n    MODEL.append(HOME / path)\n\nEVENT_MAP = {'onset':0, 'wakeup':1}\n\n# EXP = ['anglez', 'enmo', 'a_b_var', 'a_a_var', 'e_b_mean', 'e_a_mean', 'pd']\nEXP = ['anglez', 'enmo', 'a_b_var', 'e_b_mean']\nOBJ = ['state']\n\nSTEP_OF_WARMUP = 12 * 60 * 6\nSTEP_OF_RATE = 2\nSTEP_OF_STRIDE = 6\nSTEP_OF_MOVING = 12 * 60 * 6\nSTEP_OF_ROLLING = 12 * 60 * 6\n\nSEC_OF_DAY  = 60 * 60 * 24\nSEC_OF_HOUR = 60 * 60\n\nSTEP_OF_DAY   = 12 * 60 * 24\nSTEP_OF_HOUR  = 12 * 60\nSTEP_OF_30MIN = 12 * 30","metadata":{"execution":{"iopub.status.busy":"2024-01-14T07:58:57.801362Z","iopub.execute_input":"2024-01-14T07:58:57.801722Z","iopub.status.idle":"2024-01-14T07:58:57.811657Z","shell.execute_reply.started":"2024-01-14T07:58:57.801691Z","shell.execute_reply":"2024-01-14T07:58:57.810764Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# データ読込","metadata":{}},{"cell_type":"code","source":"%%time\n\ndf = pq.read_table(CUSTOM / 'custom_train_series.parquet').to_pandas()\ndf_test = pq.read_table(ORIGIN / 'test_series.parquet').to_pandas()\ndf_event = pd.read_csv(CUSTOM / 'custom_train_events.csv', index_col=0)\ndf_sample = pd.read_csv(ORIGIN / 'sample_submission.csv', index_col=0)\ndf_id = pd.read_csv(CUSTOM / 'custom_train_stat_groupby_id.csv', index_col=0)\ndf_mean = pd.read_csv(CUSTOM / 'custom_train_series_mean.csv', index_col=0)\ndf_metric = pd.read_csv(CUSTOM / 'metric.csv', index_col=0)","metadata":{"execution":{"iopub.status.busy":"2024-01-14T07:58:57.813032Z","iopub.execute_input":"2024-01-14T07:58:57.813370Z","iopub.status.idle":"2024-01-14T07:59:33.954270Z","shell.execute_reply.started":"2024-01-14T07:58:57.813326Z","shell.execute_reply":"2024-01-14T07:59:33.953250Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"series_id = df_id[['series_id']].to_dict()['series_id']","metadata":{"execution":{"iopub.status.busy":"2024-01-14T07:59:33.955301Z","iopub.execute_input":"2024-01-14T07:59:33.955595Z","iopub.status.idle":"2024-01-14T07:59:33.966295Z","shell.execute_reply.started":"2024-01-14T07:59:33.955554Z","shell.execute_reply":"2024-01-14T07:59:33.965434Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def data_pre(df):\n    df['time'] = df['timestamp'].apply(lambda x: x[:-5]).astype('datetime64[ns]')\n    df = df.merge(df_id[['series_id', 'id_map']], on='series_id', how=\"left\").drop(['series_id', 'timestamp'], axis=1)\n    df.info()\n    return df","metadata":{"execution":{"iopub.status.busy":"2024-01-14T07:59:33.967551Z","iopub.execute_input":"2024-01-14T07:59:33.968019Z","iopub.status.idle":"2024-01-14T07:59:33.977672Z","shell.execute_reply.started":"2024-01-14T07:59:33.967985Z","shell.execute_reply":"2024-01-14T07:59:33.976934Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_test = data_pre(df_test)","metadata":{"execution":{"iopub.status.busy":"2024-01-14T07:59:33.978691Z","iopub.execute_input":"2024-01-14T07:59:33.978941Z","iopub.status.idle":"2024-01-14T07:59:34.024850Z","shell.execute_reply.started":"2024-01-14T07:59:33.978919Z","shell.execute_reply":"2024-01-14T07:59:34.023948Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def moving_var(x_array, term=None, w=STEP_OF_MOVING):\n    w = min(w, len(x_array)//4)\n    m_init = np.mean(x_array[:w]); mean_list = [m_init] #平均の初期値\n    v_init = np.var(x_array[:w]); var_list = [v_init] #分散の初期値\n    for i in range(len(x_array)-w): # 更新式\n        mean_list.append(mean_list[-1] + (-x_array[i] +x_array[i+w])/w)\n        var_list.append(var_list[-1]+mean_list[-2]**2 +(x_array[i+w]**2 - x_array[i]**2)/w  -mean_list[-1]**2)\n    if term == 'before':\n        mean_list = mean_list + [mean_list[-1] for _ in range(w-1)]        \n        var_list = var_list + [var_list[-1] for _ in range(w-1)]\n    if term == 'after':\n        mean_list = [mean_list[0] for _ in range(w-1)] + mean_list\n        var_list = [var_list[0] for _ in range(w-1)] + var_list\n    return mean_list, var_list","metadata":{"execution":{"iopub.status.busy":"2024-01-14T07:59:34.029792Z","iopub.execute_input":"2024-01-14T07:59:34.030063Z","iopub.status.idle":"2024-01-14T07:59:34.039339Z","shell.execute_reply.started":"2024-01-14T07:59:34.030040Z","shell.execute_reply":"2024-01-14T07:59:34.038269Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"graphs(df, 0)","metadata":{"execution":{"iopub.status.busy":"2024-01-14T07:59:34.040362Z","iopub.execute_input":"2024-01-14T07:59:34.040677Z","iopub.status.idle":"2024-01-14T07:59:38.629861Z","shell.execute_reply.started":"2024-01-14T07:59:34.040643Z","shell.execute_reply":"2024-01-14T07:59:38.628766Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_tmp = df[(df['id_map']==0)&(df['state']!=-1)].reset_index(drop=True)\ndf_tmp","metadata":{"execution":{"iopub.status.busy":"2024-01-14T07:59:38.631341Z","iopub.execute_input":"2024-01-14T07:59:38.631747Z","iopub.status.idle":"2024-01-14T07:59:38.951097Z","shell.execute_reply.started":"2024-01-14T07:59:38.631712Z","shell.execute_reply":"2024-01-14T07:59:38.950128Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"m1, v1 = moving_var(df_tmp['anglez'], 'after', 12*60)\nm3, v3 = moving_var(df_tmp['anglez'], 'after', 12*60*3)\nm6, v6 = moving_var(df_tmp['anglez'], 'after', 12*60*6)\n\n_, axes = plt.subplots(3, 1, figsize=(20,8))\naxes[0].plot(v1)\naxes[1].plot(v3)\naxes[2].plot(v6)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-01-14T07:59:38.952256Z","iopub.execute_input":"2024-01-14T07:59:38.952559Z","iopub.status.idle":"2024-01-14T08:00:10.561221Z","shell.execute_reply.started":"2024-01-14T07:59:38.952532Z","shell.execute_reply":"2024-01-14T08:00:10.560270Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"m1, v1 = moving_var(df_tmp['enmo'], 'after', 12*60)\nm3, v3 = moving_var(df_tmp['enmo'], 'after', 12*60*3)\nm6, v6 = moving_var(df_tmp['enmo'], 'after', 12*60*6)\n\n_, axes = plt.subplots(3, 1, figsize=(20,8))\naxes[0].plot(m1)\naxes[1].plot(m3)\naxes[2].plot(m6)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:00:10.562543Z","iopub.execute_input":"2024-01-14T08:00:10.562852Z","iopub.status.idle":"2024-01-14T08:00:42.577504Z","shell.execute_reply.started":"2024-01-14T08:00:10.562825Z","shell.execute_reply":"2024-01-14T08:00:42.576623Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 学習データの加工","metadata":{}},{"cell_type":"code","source":"def add_pd(df, id_map, df_mean=df_mean):\n    df_tmp = df_mean[df_mean['id_map']==id_map][['time', 'w_pd', 'o_pd']].rename(columns={'time':'t'})\n    df['t'] = df['time'].apply(time_to_step).astype('int64')\n    df = df.merge(df_tmp, on='t').sort_values('step').drop('t', axis=1)\n    df['pd'] = df['w_pd'] - df['o_pd']\n    df = df.drop(['o_pd', 'w_pd'], axis=1)\n    return df\n\ndef add_warmup(df, id_map, warmup, col, df_mean=df_mean):\n    if warmup==0: return df\n\n    ts = df['time'].iat[0]\n    ts = ts.timestamp() % SEC_OF_DAY / 5\n    \n    df_tmp = df_mean[df_mean['id_map']==id_map]\n    df_tmp = pd.concat([df_tmp[df_tmp['time']>=ts], df_tmp[df_tmp['time']<ts]])    \n    df_tmp = df_tmp.iloc[-warmup:,:].reset_index(drop=True)\n    df_tmp['step'] = df_tmp.index - warmup\n    df_tmp['pd'] = df_tmp['w_pd'] - df_tmp['o_pd']\n    df_tmp = df_tmp.drop(['o_pd', 'w_pd'], axis=1)\n\n    ts = df['time'].iat[0]\n    ts = [ts + timedelta(seconds=-5*(i+1)) for i in reversed(range(warmup))]\n    df_tmp['time'] = ts\n    df = pd.concat([df_tmp[col], df[col]])\n    return df","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:00:42.578743Z","iopub.execute_input":"2024-01-14T08:00:42.579134Z","iopub.status.idle":"2024-01-14T08:00:42.590975Z","shell.execute_reply.started":"2024-01-14T08:00:42.579098Z","shell.execute_reply":"2024-01-14T08:00:42.590251Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def train_test_extract(id_map, warmup, n_val=2, n_test=2, df=df, df_id=df_id):\n    df = df[(df['id_map']==id_map) & (df['state']!=-1)].copy() if 'state' in df.columns else df.copy()\n    df = add_pd(df, id_map).reset_index(drop=True)\n    \n    night   = df_id[df_id['id_map']==id_map]['night'].iat[0]\n    n_train = night - n_test - n_val\n    steps = df_event[(df_event['id_map']==id_map) & (df_event['step']!=-1)]['step']\n\n    if n_test:\n        end    = steps.iat[(n_train+n_val)*2-1]\n        start  = steps.iat[(n_train+n_val)*2]\n        s_test = (end + start) // 2\n    else:\n        s_test = df['step'].iat[-1] + 1\n        \n    if n_val:\n        end   = steps.iat[n_train*2-1]\n        start = steps.iat[n_train*2]\n        s_val = (end + start) // 2\n    else:\n        s_val = s_test\n\n    col = ['id_map', 'step', 'time'] + EXP + (OBJ if 'state' in df.columns else [])\n    df = df[col]\n        \n    df_tmp = df[df['step']<s_val].reset_index(drop=True)\n    train = add_warmup(df_tmp, id_map, warmup, col)\n    df_tmp = df[(df['step']>=s_val) & (df['step']<s_test)].reset_index(drop=True)\n    val = add_warmup(df_tmp, id_map, warmup, col) if n_val else df.head(0)\n    df_tmp = df[df['step']>=s_test].reset_index(drop=True)\n    test  = add_warmup(df_tmp, id_map, warmup, col) if n_test else df.head(0)\n    \n    IPython.display.clear_output()\n    print(f'night: {night}')\n    print(f'train: {n_train}')\n    print(f'val  : {n_val}')\n    print(f'test : {n_test}')\n    \n    return train, val, test","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:00:42.592073Z","iopub.execute_input":"2024-01-14T08:00:42.592349Z","iopub.status.idle":"2024-01-14T08:00:42.613258Z","shell.execute_reply.started":"2024-01-14T08:00:42.592325Z","shell.execute_reply":"2024-01-14T08:00:42.612548Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def calc_window(hour):\n    rolling = 12 * 60 * hour // STEP_OF_STRIDE\n    warmup = STEP_OF_WARMUP + rolling * STEP_OF_STRIDE\n    return rolling, warmup","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:00:42.614233Z","iopub.execute_input":"2024-01-14T08:00:42.614496Z","iopub.status.idle":"2024-01-14T08:00:42.631137Z","shell.execute_reply.started":"2024-01-14T08:00:42.614472Z","shell.execute_reply":"2024-01-14T08:00:42.630228Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# モデルの学習","metadata":{}},{"cell_type":"code","source":"warmup = STEP_OF_WARMUP + STEP_OF_ROLLING # stride処理前のstepに対してのwarmup\nstride = STEP_OF_STRIDE\nrate   = STEP_OF_RATE\n\nprint(f'warmup: {warmup}')\nprint(f'stride: {stride}')\nprint(f'rate  : {rate}')","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:00:42.632295Z","iopub.execute_input":"2024-01-14T08:00:42.633085Z","iopub.status.idle":"2024-01-14T08:00:42.643011Z","shell.execute_reply.started":"2024-01-14T08:00:42.633052Z","shell.execute_reply":"2024-01-14T08:00:42.642247Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def std_arr(arr):\n    return (arr - arr.mean()) / arr.std()\n\ndef std_df(df, cols=EXP):\n    df = df.copy()\n    for col in cols:\n        mean = df[col].mean()\n        std = df[col].std()\n        df[col] = (df[col] - mean) / std\n    return df\n\ndef make_dataset(df, stride=STEP_OF_STRIDE, rate=STEP_OF_RATE, warmup=STEP_OF_WARMUP):\n    ds = timeseries_dataset_from_array(\n        data=df[EXP],\n        targets=df[OBJ][warmup:] if OBJ[0] in df.columns else None,\n        sequence_length=warmup // rate,\n        sampling_rate=rate,\n        sequence_stride=stride\n    )\n    return ds\n\ndef class_weight(df):\n    weight = len(df[df['state']==1]) / len(df[df['state']==0])\n    return {0: weight, 1: 1}\n\ndef handle_history(history, id_map):\n    df = pd.DataFrame.from_dict(history)\n    df['epoch'] = df.index\n    df['id_map'] = id_map\n    return df[['id_map','epoch','loss','accuracy','val_loss','val_accuracy']]\n\ndef prediction(model, df, id_map=-1, stride=STEP_OF_STRIDE, warmup=STEP_OF_WARMUP):\n    ds = make_dataset(std_df(df), stride)\n    pred = model.predict(ds).flatten()\n    pred = pd.Series(pred, index=range(0, len(pred)*stride, stride))\n    plt.figure(figsize=(20,4))\n    plt.plot(pred)\n    plt.plot(df.iloc[warmup:,:]['state'].reset_index(drop=True))\n    plt.title(f\"model: {id_map:03d}, pred: {df['id_map'].iat[0]:03d}\")\n    plt.show()\n    return pred\n\ndef predictions(models, df, stride=STEP_OF_STRIDE, warmup=STEP_OF_WARMUP):\n    preds = pd.Series()\n    pred_ = []\n    ds = make_dataset(std_df(df), stride)\n    for n_model, (id_map, model) in enumerate(models.items()):\n        pred = model.predict(ds).flatten()\n        pred_ = pred_ + pred if len(pred_) else pred\n        preds[n_model+1] = pd.Series(pred_/(n_model+1), index=range(0, len(pred)*stride, stride))\n    return preds","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:18:16.873885Z","iopub.execute_input":"2024-01-14T08:18:16.874652Z","iopub.status.idle":"2024-01-14T08:18:16.890872Z","shell.execute_reply.started":"2024-01-14T08:18:16.874617Z","shell.execute_reply":"2024-01-14T08:18:16.889592Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def handle_model(id_map, key=False ,model=False):\n    if not model:\n        model = Sequential([\n            LSTM(64, return_sequences=False),\n            Dropout(0.1),\n            Dense(units=1, activation='sigmoid')\n        ])\n    model.compile(optimizer=Adam(),\n        loss=BinaryCrossentropy(from_logits=False),\n        metrics=[\"binary_crossentropy\"])\n    \n    early_stopping = EarlyStopping(monitor='val_loss', patience=5, mode='min')\n\n    key = ('_'+key) if key else ''\n    checkpoint_callback = ModelCheckpoint(\n        filepath=f'/kaggle/working/model_{id_map:03d}{key}.tf',\n        monitor='val_loss',\n        mode='min',\n        save_best_only=True,\n        save_weights_only=False,\n        verbose=1\n    )\n    \n    return model, early_stopping, checkpoint_callback","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:00:42.660874Z","iopub.execute_input":"2024-01-14T08:00:42.661192Z","iopub.status.idle":"2024-01-14T08:00:42.683314Z","shell.execute_reply.started":"2024-01-14T08:00:42.661162Z","shell.execute_reply":"2024-01-14T08:00:42.682620Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(f\"モデル作成用: {len(df_id[df_id['night']>=20]['night'])}\")\nprint(f\"検証用　　　: {len(df_id[(df_id['night']<20)&(df_id['night']!=0)]['night'])}\")","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:00:42.684321Z","iopub.execute_input":"2024-01-14T08:00:42.684656Z","iopub.status.idle":"2024-01-14T08:00:42.703262Z","shell.execute_reply.started":"2024-01-14T08:00:42.684629Z","shell.execute_reply":"2024-01-14T08:00:42.702424Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"models = pd.Series()","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:00:42.704625Z","iopub.execute_input":"2024-01-14T08:00:42.704987Z","iopub.status.idle":"2024-01-14T08:00:42.713620Z","shell.execute_reply.started":"2024-01-14T08:00:42.704953Z","shell.execute_reply":"2024-01-14T08:00:42.712945Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"c = 0\nfor i, id_map in enumerate(df_id[df_id['night']>=20].index):\n    if i % 50 == 0:\n        c += 1\n    path = f\"input/zzz-model-{c}-ipynb/model_{id_map:03d}.tf\"\n    models[id_map] = tf.keras.models.load_model(HOME / path)","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:00:42.714753Z","iopub.execute_input":"2024-01-14T08:00:42.715116Z","iopub.status.idle":"2024-01-14T08:05:11.978706Z","shell.execute_reply.started":"2024-01-14T08:00:42.715082Z","shell.execute_reply":"2024-01-14T08:05:11.977797Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"models","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:05:11.979845Z","iopub.execute_input":"2024-01-14T08:05:11.980108Z","iopub.status.idle":"2024-01-14T08:05:11.988115Z","shell.execute_reply.started":"2024-01-14T08:05:11.980085Z","shell.execute_reply":"2024-01-14T08:05:11.987193Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_pre_score = []\nfor i in range(1, 4):\n    path = f\"input/zzz-model-{i}-ipynb/pred_score.csv\"\n    df_tmp = pd.read_csv(HOME / path, index_col=0)\n    df_pre_score = pd.concat([df_pre_score, df_tmp]) if len(df_pre_score) else df_tmp\ndf_pre_score = df_pre_score.reset_index(drop=True)\ndf_pre_score.sort_values('F1', ascending=False).head(10)","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:05:11.989197Z","iopub.execute_input":"2024-01-14T08:05:11.989468Z","iopub.status.idle":"2024-01-14T08:05:12.053430Z","shell.execute_reply.started":"2024-01-14T08:05:11.989445Z","shell.execute_reply":"2024-01-14T08:05:12.052493Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.hist(df_pre_score['F1'])\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:05:12.054727Z","iopub.execute_input":"2024-01-14T08:05:12.055047Z","iopub.status.idle":"2024-01-14T08:05:12.321579Z","shell.execute_reply.started":"2024-01-14T08:05:12.055020Z","shell.execute_reply":"2024-01-14T08:05:12.320520Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 変化点の検出","metadata":{}},{"cell_type":"code","source":"def pred_pre(pred, rolling=STEP_OF_ROLLING, stride=STEP_OF_STRIDE):\n    w = STEP_OF_ROLLING // STEP_OF_STRIDE\n    pred = pred.rolling(w).mean()[w-1:].reset_index(drop=True)\n    pred = [[i*stride, i, x-0.5] for i, x in enumerate(pred)]\n    pred = pd.DataFrame(pred, columns=['idx', 'sequence', 'state'])\n    return pred","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:05:12.323068Z","iopub.execute_input":"2024-01-14T08:05:12.323947Z","iopub.status.idle":"2024-01-14T08:05:12.330351Z","shell.execute_reply.started":"2024-01-14T08:05:12.323906Z","shell.execute_reply":"2024-01-14T08:05:12.329391Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def event_handle(pred, test, stride=STEP_OF_STRIDE, series_id=series_id):\n    series_id = series_id[test['id_map'].iat[0]]\n    d = 12 * 30\n    w = 12 * 60 * 18\n    step = 0\n    state = 'onset'\n    event = []\n    while True:\n        s = min(step + d, pred['idx'].iat[-1])\n        e = min(s + w, pred['idx'].iat[-1])\n        p_tmp = pred[(pred['idx']>=s) & (pred['idx']<=e)]\n        if state == 'onset':\n            i = p_tmp['state'].idxmax()\n        else:\n            i = p_tmp['state'].idxmin()\n        step = p_tmp[p_tmp.index==i]['idx'].iat[0]\n        event.append([series_id, step, state, 0.5])\n        state = 'onset' if state=='wakeup' else 'wakeup'\n        if (state=='onset') & (e==pred['idx'].iat[-1]):\n            break\n    return event","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:05:12.335578Z","iopub.execute_input":"2024-01-14T08:05:12.335963Z","iopub.status.idle":"2024-01-14T08:05:12.345512Z","shell.execute_reply.started":"2024-01-14T08:05:12.335923Z","shell.execute_reply":"2024-01-14T08:05:12.344616Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def check_event(pred, event):\n    plt.figure(figsize=(10,4))\n    plt.plot(p_tmp['idx'], p_tmp['state'])\n    for x in event:\n        plt.axvline(x[1], c='r', linewidth=0.5)","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:05:12.346554Z","iopub.execute_input":"2024-01-14T08:05:12.346903Z","iopub.status.idle":"2024-01-14T08:05:12.371694Z","shell.execute_reply.started":"2024-01-14T08:05:12.346874Z","shell.execute_reply":"2024-01-14T08:05:12.370536Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def state_handle(df, event):\n    df = df[df['step']>=0].copy()\n    state = [1 for _ in range(event[0][1])]\n    for i in range(1, len(event)):\n        x = 1 if event[i][2] == 'onset' else 0\n        state += [x for _ in range(event[i][1]-event[i-1][1])]\n    state += [1 for _ in range(event[-1][1], len(df))]\n    df['pred'] = state\n    df.index.name = ''\n    return df","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:05:12.373171Z","iopub.execute_input":"2024-01-14T08:05:12.373582Z","iopub.status.idle":"2024-01-14T08:05:12.383773Z","shell.execute_reply.started":"2024-01-14T08:05:12.373531Z","shell.execute_reply":"2024-01-14T08:05:12.382637Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def prob_handle(df, pred, stride=STEP_OF_STRIDE):\n    prob = np.array([[x for _ in range(stride)] for x in pred['state']]).flatten()[-len(df):]\n    df['prob'] = prob\n    return df","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:05:12.385390Z","iopub.execute_input":"2024-01-14T08:05:12.385725Z","iopub.status.idle":"2024-01-14T08:05:12.396848Z","shell.execute_reply.started":"2024-01-14T08:05:12.385696Z","shell.execute_reply":"2024-01-14T08:05:12.395780Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"id_ = df_pre_score.sort_values('F1', ascending=False).head(6)['model'].unique()\nmodels_ = pd.Series()\nfor i in id_:\n    models_[i] = models[models.index==i].loc[i]\nmodels_","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:30:40.656498Z","iopub.execute_input":"2024-01-14T08:30:40.657414Z","iopub.status.idle":"2024-01-14T08:30:40.672140Z","shell.execute_reply.started":"2024-01-14T08:30:40.657376Z","shell.execute_reply":"2024-01-14T08:30:40.671098Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"columns = ['model','target','Accuracy','Recall','Precision','F1','TP','FN','FP','TN','n_TP','n_FN','n_FP','n_TN']\ndf_score = pd.DataFrame(columns=columns)\ndf_score","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:29:52.744832Z","iopub.execute_input":"2024-01-14T08:29:52.745213Z","iopub.status.idle":"2024-01-14T08:29:52.759482Z","shell.execute_reply.started":"2024-01-14T08:29:52.745184Z","shell.execute_reply":"2024-01-14T08:29:52.758415Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for target in tqdm(df_id[(df_id['night']<20)&(df_id['night']!=0)]['id_map'].unique()):\n    test, _, _ = train_test_extract(target, warmup, 0, 0)\n    preds = predictions(models_, test)\n    for n_model, pred in preds.items():\n        pred = pred_pre(pred)\n        event = event_handle(pred, test[test['step']>=0])\n        result = state_handle(test, event)\n        result = prob_handle(result, pred)\n\n        accuracy = accuracy_score(result['state'], result['pred'])\n        recall = recall_score(result['state'], result['pred'])\n        precision = precision_score(result['state'], result['pred'])\n        f1 = f1_score(result['state'], result['pred'])\n        cmx_data = confusion_matrix(result['state'], result['pred']).flatten()\n\n        score = [n_model, test['id_map'].iat[0]]\n        score += [accuracy, recall, precision, f1]\n        score += list(cmx_data / np.sum(cmx_data))\n        score += list(cmx_data)\n        df_score.loc[len(df_score)] = score\n        result[['step','state','pred','prob']].to_csv(f'pred_{target:03d}_{n_model}.csv')\n\nfor col in ['model','target','n_TP','n_FN','n_FP','n_TN']:\n    df_score[col] = df_score[col].astype(int)","metadata":{"execution":{"iopub.status.busy":"2024-01-14T08:30:58.955303Z","iopub.execute_input":"2024-01-14T08:30:58.956251Z","iopub.status.idle":"2024-01-14T08:32:25.579525Z","shell.execute_reply.started":"2024-01-14T08:30:58.956214Z","shell.execute_reply":"2024-01-14T08:32:25.578318Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_score.to_csv('pred_score.csv')\ndf_score","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# モデルの評価","metadata":{}},{"cell_type":"code","source":"# 実測-予測\n# 睡眠-睡眠：True Positive (TP)\n# 睡眠-起床：False Negative (FN)\n# 起床-睡眠：False Positive (FP)\n# 起床-起床：True Negative (TN)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_score.info()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"d = {\n    'Accuracy' : 'Accuracy（正解率）\\nどれだけ正確に予測できているか',\n    'Recall'   : 'Recall（再現率）\\nどれだけ取りこぼしなく予測することができたか',\n    'Precision': 'Precision（適合率）\\n正と予測したものがどれだけ正しかったか',\n    'F1'       : 'F1\\nRecallとPrecisionのバランス'}\n_, axes = plt.subplots(2, 2, figsize=(20,10))\naxes = axes.ravel()\nfor i, (k, v) in enumerate(d.items()):\n    axes[i].hist(df_score[k])\n    axes[i].set_title(v)\nplt.show()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"d = {\n    'TP': 'True Positive (TP, 真陽性)\\n睡眠(実)-睡眠(予)',\n    'FN': 'False Negative (FN, 偽陰性)\\n睡眠(実)-起床(予)',\n    'FP': 'False Positive (FP, 偽陽性)\\n起床(実)-睡眠(予)',\n    'TN': 'True Negative (TN, 真陰性)\\n起床(実)-起床(予)'}\n_, axes = plt.subplots(2, 2, figsize=(20,10))\naxes = axes.ravel()\nfor i, (k, v) in enumerate(d.items()):\n    axes[i].hist(df_score[k])\n    axes[i].set_title(v)\nplt.show()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"cmx_data = confusion_matrix(result['state'], result['pred'])\ncmx_data = cmx_data / np.sum(cmx_data)\nsns.heatmap(cmx_data, annot=True, cmap='Blues', fmt='.2f')\nplt.xlabel('予測')\nplt.ylabel('実測')\nplt.xticks([0.5,1.5], ['睡眠','起床'])\nplt.yticks([0.5,1.5], ['睡眠','起床'])\nplt.show()","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}