{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.14","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"nvidiaTeslaT4","dataSources":[{"sourceId":84493,"databundleVersionId":9871156,"sourceType":"competition"}],"dockerImageVersionId":30786,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport polars as pl\nfrom matplotlib import pyplot as plt\nimport seaborn as sns\nfrom pathlib import Path\nfrom tqdm import tqdm\nfrom functools import reduce\n\nfrom scipy.stats import pearsonr, spearmanr, kendalltau\n\nDATA_DIR = Path('/kaggle/input/jane-street-real-time-market-data-forecasting')\nN_PARTITION = 10\n\nfeature_cols = [f'feature_{x:02}' for x in range(79)]\nresponder_cols = [f'responder_{i}' for i in range(9)]","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2024-12-27T01:07:54.441470Z","iopub.execute_input":"2024-12-27T01:07:54.441732Z","iopub.status.idle":"2024-12-27T01:07:57.228855Z","shell.execute_reply.started":"2024-12-27T01:07:54.441705Z","shell.execute_reply":"2024-12-27T01:07:57.228000Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Prepare data and have some basic EDA","metadata":{}},{"cell_type":"code","source":"train_parquets = [DATA_DIR / f\"train.parquet/partition_id={i}/part-0.parquet\" for i in range(N_PARTITION)]\n\npartition_dates = {}\nfor _f in train_parquets:\n    partition_id = int(_f.parents[0].stem.split('=')[-1])\n    _pl = pl.read_parquet(_f, columns=[\"date_id\", \"time_id\"])\n    partition_dates[partition_id] = (_pl['date_id'].min(), _pl['date_id'].max())\n\n# check date range in each partition\npartition_dates","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-27T01:07:57.234027Z","iopub.execute_input":"2024-12-27T01:07:57.234348Z","iopub.status.idle":"2024-12-27T01:07:57.753810Z","shell.execute_reply.started":"2024-12-27T01:07:57.234310Z","shell.execute_reply":"2024-12-27T01:07:57.752914Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"pl_train = pl.read_parquet(train_parquets)\npl_train.head()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-27T01:07:59.912392Z","iopub.execute_input":"2024-12-27T01:07:59.913187Z","iopub.status.idle":"2024-12-27T01:08:43.095420Z","shell.execute_reply.started":"2024-12-27T01:07:59.913138Z","shell.execute_reply":"2024-12-27T01:08:43.094563Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"pl_train.tail()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T06:32:18.741188Z","iopub.execute_input":"2024-12-26T06:32:18.741691Z","iopub.status.idle":"2024-12-26T06:32:18.756302Z","shell.execute_reply.started":"2024-12-26T06:32:18.741648Z","shell.execute_reply":"2024-12-26T06:32:18.754862Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import gc\nimport torch\nimport torch.nn as nn\nfrom torch.utils.data import DataLoader, TensorDataset\nimport numpy as np\nimport matplotlib.pyplot as plt\nfrom sklearn.preprocessing import MinMaxScaler\n\n# LSTM 모델 정의\nclass LSTMModel(nn.Module):\n    def __init__(self, input_size, hidden_size, output_size, num_layers=1):\n        super(LSTMModel, self).__init__()\n        self.lstm = nn.LSTM(input_size, hidden_size, num_layers, batch_first=True)\n        self.fc = nn.Linear(hidden_size, output_size)\n\n    def forward(self, x):\n        _, (hn, _) = self.lstm(x)  # 마지막 hidden state 사용\n        out = self.fc(hn[-1])\n        return out\n\n# 가중치 초기화 함수\ndef init_weights(m):\n    if isinstance(m, nn.LSTM):\n        for name, param in m.named_parameters():\n            if 'weight' in name:\n                nn.init.xavier_uniform_(param)\n            elif 'bias' in name:\n                nn.init.zeros_(param)\n\n# 데이터 전처리 함수\ndef preprocess_chunk(chunk, features, target, scaler):\n    chunk = chunk.fill_nan(0)  # 결측값 처리\n    X = chunk[features].to_numpy()\n    y = chunk[target].to_numpy()\n    X_scaled = scaler.transform(X)  # 스케일링\n    X_tensor = torch.tensor(X_scaled, dtype=torch.float32).unsqueeze(1).to(device)\n    y_tensor = torch.tensor(y, dtype=torch.float32).to(device)\n\n    # NaN 또는 Inf 제거\n    X_tensor = torch.nan_to_num(X_tensor, nan=0.0, posinf=1e6, neginf=-1e6)\n    y_tensor = torch.nan_to_num(y_tensor, nan=0.0, posinf=1e6, neginf=-1e6)\n    \n    return X_tensor, y_tensor\n\n# 학습 함수\ndef train_model_on_chunks(pl_train, features, target, chunk_size, model, optimizer, criterion, scaler):\n    model.train()\n    for start in range(0, len(pl_train), chunk_size):\n        end = start + chunk_size\n        chunk = pl_train[start:end]\n        X_chunk, y_chunk = preprocess_chunk(chunk, features, target, scaler)\n\n        optimizer.zero_grad()\n        output = model(X_chunk)\n        loss = criterion(output.squeeze(), y_chunk)\n\n        if torch.isnan(loss):\n            print(\"NaN detected in loss\")\n            print(\"Output:\", output[:5].detach().cpu().numpy())\n            print(\"Target:\", y_chunk[:5].cpu().numpy())\n            break\n\n        loss.backward()\n        torch.nn.utils.clip_grad_norm_(model.parameters(), max_norm=1.0)  # Gradient Clipping\n        optimizer.step()\n\n        print(f\"Chunk {start // chunk_size + 1}: Loss = {loss.item()}\")\n        del chunk, X_chunk, y_chunk\n        gc.collect()\n        torch.cuda.empty_cache()\n\n# GPU 초기화\ndevice = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\n\n# 데이터 준비\nfeatures = [f\"feature_{i:02d}\" for i in range(79)]\ntarget = \"responder_6\"\nscaler = MinMaxScaler()\n\n# 모델 초기화\ninput_size = len(features)\nhidden_size = 64\noutput_size = 1\nnum_layers = 2\nmodel = LSTMModel(input_size, hidden_size, output_size, num_layers).to(device)\nmodel.apply(init_weights)\noptimizer = torch.optim.Adam(model.parameters(), lr=0.001)\ncriterion = nn.MSELoss()\n\n# 스케일러 학습\nsample_data = pl_train[:5000]  # 일부 데이터로 스케일러 학습\nscaler.fit(sample_data[features].to_numpy())\n\n# 학습 실행\nchunk_size = 10000  # 청크 크기\ntrain_model_on_chunks(pl_train, features, target, chunk_size, model, optimizer, criterion, scaler)\n\n# 예측\nmodel.eval()\nwith torch.no_grad():\n    sample_chunk = pl_train[:5000]  # 일부 데이터로 테스트\n    X_test, y_test = preprocess_chunk(sample_chunk, features, target, scaler)\n    predictions = model(X_test).cpu().numpy()\n\n# 결과 시각화\nplt.figure(figsize=(12, 6))\nplt.plot(y_test.cpu().numpy(), label=\"True Values\", alpha=0.7)\nplt.plot(predictions, label=\"Predictions\", alpha=0.7)\nplt.legend()\nplt.show()\n\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-27T01:45:19.967801Z","iopub.execute_input":"2024-12-27T01:45:19.968145Z","iopub.status.idle":"2024-12-27T01:59:39.647840Z","shell.execute_reply.started":"2024-12-27T01:45:19.968115Z","shell.execute_reply":"2024-12-27T01:59:39.647022Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# 해석\n초기 결측치:\n\n초기 date_id 구간에서는 일부 feature에서 결측치가 집중적으로 발생합니다.\n특히 feature_00 ~ feature_04는 초반에 결측치 비율이 매우 높습니다.\n후반부 결측치 감소:\n\n시간이 지남에 따라 대부분의 feature에서 결측치가 사라지거나 감소하는 경향이 보입니다.\n이는 데이터 수집 시스템이 점진적으로 개선되었음을 의미합니다.\n특이 feature:\n\n일부 feature는 초반에 존재하지 않다가 중간 또는 후반부에 데이터가 추가된 것을 확인할 수 있습니다.\n\n# 분석 활용\n결측치 처리:\n결측치가 많은 feature는 분석 대상에서 제외하거나 보간법 등을 사용해 처리할 수 있습니다.\n데이터 수집 개선:\n결측치 발생 시점 및 패턴을 통해 데이터 수집 과정에서 발생한 문제점을 파악합니다.\n시계열 모델링:\n결측치가 존재하는 기간과 feature를 구분하여 모델링 시 결측치 보정이 필요합니다.","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nimport matplotlib.pyplot as plt\n\n# 데이터 로드 및 분할\npl_train = pl.read_parquet(train_parquets, columns=[\"date_id\", \"time_id\", \"responder_6\"])\n\n# date_id 676을 기준으로 데이터 구분\nearly_data = pl_train.filter(pl.col(\"date_id\") <= 676)\nlate_data = pl_train.filter(pl.col(\"date_id\") > 676)\n\n# 보간 함수 (중복 time_id 처리 추가)\ndef interpolate_data(pl_data, time_range):\n    df = pl_data.to_pandas()\n    df = df.sort_values(['date_id', 'time_id'])\n    \n    interpolated_data = []\n    for date in df['date_id'].unique():\n        subset = df[df['date_id'] == date]\n        \n        # 중복된 time_id 처리 (평균값을 사용)\n        subset = subset.groupby('time_id')['responder_6'].mean().reset_index()\n        subset = subset.set_index('time_id')\n\n        # 새로운 시간 인덱스 생성 (리샘플링)\n        new_index = pd.RangeIndex(start=time_range[0], stop=time_range[1] + 1)\n        subset = subset.reindex(new_index)\n\n        # 선형 보간 수행\n        subset = subset.interpolate(method='linear')\n        interpolated_data.append(subset.reset_index().assign(date_id=date))\n    \n    return pd.concat(interpolated_data)\n\n# 676 이전 -> 0~967로 리샘플링\nearly_interpolated = interpolate_data(early_data, (0, 967))\n# 676 이후 -> 그대로 사용\nlate_interpolated = interpolate_data(late_data, (0, 967))\n\n# 병합\nfinal_interpolated = pd.concat([early_interpolated, late_interpolated])\nfinal_interpolated = final_interpolated.rename(columns={\"index\": \"time_id\"})\n\n# 시각화\ndate_choice = 500  # 676 이전\ndate_choice_late = 700  # 676 이후\n\nfig, axes = plt.subplots(1, 2, figsize=(14, 6), sharey=True)\n\n# 676 이전 데이터 시각화\nsubset_early = final_interpolated[final_interpolated['date_id'] == date_choice]\naxes[0].plot(subset_early['time_id'], subset_early['responder_6'], label=\"Interpolated Data\", linestyle='-', alpha=0.8)\naxes[0].set_title(f\"Interpolated Responder_6 (Date ID {date_choice})\")\naxes[0].grid(True, linestyle=\"--\")\n\n# 676 이후 데이터 시각화\nsubset_late = final_interpolated[final_interpolated['date_id'] == date_choice_late]\naxes[1].plot(subset_late['time_id'], subset_late['responder_6'], label=\"Interpolated Data\", linestyle='-', alpha=0.8)\naxes[1].set_title(f\"Interpolated Responder_6 (Date ID {date_choice_late})\")\naxes[1].grid(True, linestyle=\"--\")\n\nplt.suptitle(\"Linear Interpolation for Different Time Ranges\")\nplt.show()\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T05:46:46.283733Z","iopub.execute_input":"2024-12-26T05:46:46.284107Z","iopub.status.idle":"2024-12-26T05:47:33.553015Z","shell.execute_reply.started":"2024-12-26T05:46:46.284072Z","shell.execute_reply":"2024-12-26T05:47:33.551569Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"타임스탬프(시간간격)\n0~676: 849 / 677~1698: 968\n간격이 다름으로 시계열 분석 위해 선형보간 작업.\n위는 보간 통한 리샘플링으로 타임스탬프 간격 일정해짐\n하지만 왼쪽 그래프 우측 보면 800 이후부터 값이 일정하게 나옴\n이유는 간격 차이로 인해 데이터가 적어서 그럼\n이를 해결하기 위해 외삽 추가(주변값)","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nimport matplotlib.pyplot as plt\n\n# 데이터 로드 및 분할\npl_train = pl.read_parquet(train_parquets, columns=[\"date_id\", \"time_id\", \"responder_6\"])\n\n# date_id 676을 기준으로 데이터 구분\nearly_data = pl_train.filter(pl.col(\"date_id\") <= 676)\nlate_data = pl_train.filter(pl.col(\"date_id\") > 676)\n\n# 보간 함수 (중복 time_id 처리 + 외삽 추가)\ndef interpolate_data(pl_data, time_range):\n    df = pl_data.to_pandas()\n    df = df.sort_values(['date_id', 'time_id'])\n    \n    interpolated_data = []\n    for date in df['date_id'].unique():\n        subset = df[df['date_id'] == date]\n        \n        # 중복된 time_id 처리 (평균값 사용)\n        subset = subset.groupby('time_id')['responder_6'].mean().reset_index()\n        subset = subset.set_index('time_id')\n\n        # 새로운 시간 인덱스 생성 (리샘플링)\n        new_index = pd.RangeIndex(start=time_range[0], stop=time_range[1] + 1)\n        subset = subset.reindex(new_index)\n\n        # 선형 보간 수행 + 외삽\n        subset = subset.interpolate(method='polynomial', order=3, limit_direction='both')\n        interpolated_data.append(subset.reset_index().assign(date_id=date))\n    \n    return pd.concat(interpolated_data)\n\n# 676 이전 -> 0~967로 리샘플링 및 보간\nearly_interpolated = interpolate_data(early_data, (0, 967))\n# 676 이후 -> 그대로 사용\nlate_interpolated = interpolate_data(late_data, (0, 967))\n\n# 병합\nfinal_interpolated = pd.concat([early_interpolated, late_interpolated])\nfinal_interpolated = final_interpolated.rename(columns={\"index\": \"time_id\"})\n\n# 시각화\ndate_choice = 500  # 676 이전\ndate_choice_late = 700  # 676 이후\n\nfig, axes = plt.subplots(1, 2, figsize=(14, 6), sharey=True)\n\n# 676 이전 데이터 시각화\nsubset_early = final_interpolated[final_interpolated['date_id'] == date_choice]\naxes[0].plot(subset_early['time_id'], subset_early['responder_6'], label=\"Interpolated Data\", linestyle='-', alpha=0.8)\naxes[0].set_title(f\"Interpolated Responder_6 (Date ID {date_choice})\")\naxes[0].grid(True, linestyle=\"--\")\n\n# 676 이후 데이터 시각화\nsubset_late = final_interpolated[final_interpolated['date_id'] == date_choice_late]\naxes[1].plot(subset_late['time_id'], subset_late['responder_6'], label=\"Interpolated Data\", linestyle='-', alpha=0.8)\naxes[1].set_title(f\"Interpolated Responder_6 (Date ID {date_choice_late})\")\naxes[1].grid(True, linestyle=\"--\")\n\nplt.suptitle(\"Improved Interpolation with Extrapolation\")\nplt.show()\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T05:47:33.554933Z","iopub.execute_input":"2024-12-26T05:47:33.555449Z","iopub.status.idle":"2024-12-26T05:48:20.399462Z","shell.execute_reply.started":"2024-12-26T05:47:33.555395Z","shell.execute_reply":"2024-12-26T05:48:20.398194Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"'''\nfrom statsmodels.tsa.stattools import adfuller\nimport matplotlib.pyplot as plt\n\n# ADF 테스트 함수 정의\ndef test_stationarity(series, title=\"Overall Time Series\"):\n    \"\"\"\n    주어진 시계열 데이터에 대해 Augmented Dickey-Fuller (ADF) 테스트를 수행.\n    \"\"\"\n    print(f\"== {title} ==\")\n    adf_result = adfuller(series, autolag=\"AIC\")\n    print(f\"ADF Statistic: {adf_result[0]}\")\n    print(f\"p-value: {adf_result[1]}\")\n    print(\"Critical Values:\")\n    for key, value in adf_result[4].items():\n        print(f\"   {key}: {value}\")\n    if adf_result[1] <= 0.05:\n        print(\"-> Series is stationary.\")\n    else:\n        print(\"-> Series is not stationary.\")\n    print(\"-\" * 50)\n\n# 전체 데이터에서 responder_6 병합\noverall_series = final_interpolated.sort_values([\"date_id\", \"time_id\"])[\"responder_6\"]\n\n# 정상성 테스트 수행\ntest_stationarity(overall_series, title=\"Responder_6 (Full Dataset)\")\n\n# 시각화\nplt.figure(figsize=(12, 6))\nplt.plot(overall_series.reset_index(drop=True), label=\"Responder_6\")\nplt.title(\"Responder_6 Time Series (Full Dataset)\")\nplt.xlabel(\"Time Index\")\nplt.ylabel(\"Responder_6\")\nplt.grid(ls=\"--\")\nplt.legend()\nplt.show()\n'''","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T05:48:20.401034Z","iopub.execute_input":"2024-12-26T05:48:20.401444Z","iopub.status.idle":"2024-12-26T05:48:20.410482Z","shell.execute_reply.started":"2024-12-26T05:48:20.401405Z","shell.execute_reply":"2024-12-26T05:48:20.409172Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"결측치 많아서 분석 안됨 -> 새로운 보간 작업 필요","metadata":{}},{"cell_type":"code","source":"#final_interpolated.isna().sum()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T05:48:20.414347Z","iopub.execute_input":"2024-12-26T05:48:20.414765Z","iopub.status.idle":"2024-12-26T05:48:20.429162Z","shell.execute_reply.started":"2024-12-26T05:48:20.414728Z","shell.execute_reply":"2024-12-26T05:48:20.427801Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"'''\nimport pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\n\n# 데이터 로드\npl_train = pl.read_parquet(train_parquets, columns=[\"date_id\", \"time_id\", \"responder_6\"])\n\n# date_id 676 기준으로 데이터 분리\nearly_data = pl_train.filter(pl.col(\"date_id\") <= 676)\nlate_data = pl_train.filter(pl.col(\"date_id\") > 676)\n\n# 보간 함수: 타임스탬프 통일 및 값 보간\ndef unify_time_intervals(pl_data, target_timestamps):\n    df = pl_data.to_pandas()\n    df = df.sort_values([\"date_id\", \"time_id\"])\n    \n    unified_data = []\n    for date in df[\"date_id\"].unique():\n        subset = df[df[\"date_id\"] == date]\n        \n        # 중복된 time_id 처리\n        subset = subset.groupby(\"time_id\")[\"responder_6\"].mean().reset_index()\n        \n        # 기존 time_id를 새로운 time_id로 매핑\n        new_time_id = np.linspace(\n            subset[\"time_id\"].min(), subset[\"time_id\"].max(), target_timestamps\n        )\n        \n        # 보간 수행\n        interpolated_responder = np.interp(new_time_id, subset[\"time_id\"], subset[\"responder_6\"])\n        \n        # 새로운 DataFrame 생성\n        unified_subset = pd.DataFrame({\n            \"time_id\": new_time_id,\n            \"responder_6\": interpolated_responder,\n            \"date_id\": date\n        })\n        unified_data.append(unified_subset)\n    \n    return pd.concat(unified_data, ignore_index=True)\n\n# 676 이전 -> 968 타임스탬프 기준으로 보간\nearly_interpolated = unify_time_intervals(early_data, late_timestamps)\n# 676 이후 -> 그대로 사용 (이미 968 타임스탬프)\nlate_interpolated = unify_time_intervals(late_data, late_timestamps)\n\n# 데이터 병합\nfinal_interpolated = pd.concat([early_interpolated, late_interpolated], ignore_index=True)\n\n# 시각화\ndate_choice_early = 500  # 676 이전\ndate_choice_late = 700  # 676 이후\n\nfig, axes = plt.subplots(1, 2, figsize=(14, 6), sharey=True)\n\n# 676 이전 데이터 시각화\nsubset_early = final_interpolated[final_interpolated[\"date_id\"] == date_choice_early]\naxes[0].plot(\n    subset_early[\"time_id\"],\n    subset_early[\"responder_6\"],\n    label=\"Unified Data\",\n    linestyle=\"-\",\n    alpha=0.8,\n)\naxes[0].set_title(f\"Unified Time Intervals: Responder_6 (Date ID {date_choice_early})\")\naxes[0].grid(True, linestyle=\"--\")\n\n# 676 이후 데이터 시각화\nsubset_late = final_interpolated[final_interpolated[\"date_id\"] == date_choice_late]\naxes[1].plot(\n    subset_late[\"time_id\"],\n    subset_late[\"responder_6\"],\n    label=\"Unified Data\",\n    linestyle=\"-\",\n    alpha=0.8,\n)\naxes[1].set_title(f\"Unified Time Intervals: Responder_6 (Date ID {date_choice_late})\")\naxes[1].grid(True, linestyle=\"--\")\n\nplt.suptitle(\"Unified Time Intervals Across All Dates\")\nplt.show()\n'''","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T05:48:20.430921Z","iopub.execute_input":"2024-12-26T05:48:20.431331Z","iopub.status.idle":"2024-12-26T05:48:20.446553Z","shell.execute_reply.started":"2024-12-26T05:48:20.431290Z","shell.execute_reply":"2024-12-26T05:48:20.445111Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#final_interpolated.isna().sum()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T05:48:20.448481Z","iopub.execute_input":"2024-12-26T05:48:20.449360Z","iopub.status.idle":"2024-12-26T05:48:20.466010Z","shell.execute_reply.started":"2024-12-26T05:48:20.449300Z","shell.execute_reply":"2024-12-26T05:48:20.464686Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#final_interpolated.shape","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T05:48:20.467624Z","iopub.execute_input":"2024-12-26T05:48:20.468183Z","iopub.status.idle":"2024-12-26T05:48:20.480995Z","shell.execute_reply.started":"2024-12-26T05:48:20.468123Z","shell.execute_reply":"2024-12-26T05:48:20.479380Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"'''\nfrom statsmodels.tsa.stattools import adfuller\nimport numpy as np\nimport pandas as pd\n\n# ADF 테스트 함수 정의\ndef test_stationarity_chunkwise(chunk, title=\"Time Series Chunk\"):\n    \"\"\"\n    청크 단위로 ADF 테스트를 수행하여 시계열 데이터의 정상성을 확인.\n    \"\"\"\n    print(f\"== {title} ==\")\n    chunk = chunk[~np.isnan(chunk)]  # 결측치 제거\n    adf_result = adfuller(chunk, autolag=\"AIC\")\n    print(f\"ADF Statistic: {adf_result[0]:.4f}\")\n    print(f\"p-value: {adf_result[1]:.4f}\")\n    print(\"Critical Values:\")\n    for key, value in adf_result[4].items():\n        print(f\"    {key}: {value:.4f}\")\n    if adf_result[1] <= 0.05:\n        print(\"\\n[Result] The series is stationary (p-value <= 0.05).\\n\")\n    else:\n        print(\"\\n[Result] The series is not stationary (p-value > 0.05).\\n\")\n\n# 청크 단위 처리 함수\ndef process_chunks(data, chunk_size):\n    \"\"\"\n    데이터프레임을 청크 단위로 분할하고 정상성 분석 준비.\n    \"\"\"\n    num_rows = data.shape[0]\n    for i in range(0, num_rows, chunk_size):\n        chunk = data.iloc[i:i + chunk_size]\n        yield chunk[\"responder_6\"].values\n\n# 청크 단위로 데이터 처리\nchunk_size = 100000  # 청크 크기 설정 (100,000개씩 처리)\nchunks = process_chunks(final_interpolated, chunk_size)\n\n# 정상성 분석 수행\nfor idx, chunk in enumerate(chunks):\n    print(f\"\\n--- Chunk {idx + 1} ---\")\n    test_stationarity_chunkwise(chunk, title=f\"Chunk {idx + 1}\")\n'''\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T05:48:20.482425Z","iopub.execute_input":"2024-12-26T05:48:20.482788Z","iopub.status.idle":"2024-12-26T05:48:20.497064Z","shell.execute_reply.started":"2024-12-26T05:48:20.482753Z","shell.execute_reply":"2024-12-26T05:48:20.495851Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"아래는 가중평균 고려한 정상성 분석","metadata":{}},{"cell_type":"code","source":"import polars as pl\nimport pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\n\n# Step 1: 데이터 로드 및 Pandas 변환\npl_train = pl.read_parquet(train_parquets, columns=[\"date_id\", \"time_id\", \"weight\", \"responder_6\"])\ndf = pl_train.to_pandas()\n\n# Step 2: Pandas에서 `time_id`별 가중평균 계산\ndf[\"weighted_sum\"] = df[\"responder_6\"] * df[\"weight\"]\nweighted_mean_df = (\n    df.groupby([\"date_id\", \"time_id\"])\n    .agg(\n        weighted_mean_responder_6=(\"weighted_sum\", \"sum\"),\n        total_weight=(\"weight\", \"sum\")\n    )\n    .reset_index()\n)\nweighted_mean_df[\"weighted_mean_responder_6\"] /= weighted_mean_df[\"total_weight\"]\n\n# Step 3: 타임스탬프 간격 맞추기\n# date_id 676 기준으로 데이터 분리\nearly_data = weighted_mean_df[weighted_mean_df[\"date_id\"] <= 676]\nlate_data = weighted_mean_df[weighted_mean_df[\"date_id\"] > 676]\n\n# 보간 및 외삽 함수\ndef interpolate_and_extrapolate(df, target_timestamps):\n    interpolated_data = []\n    for date in df[\"date_id\"].unique():\n        subset = df[df[\"date_id\"] == date]\n        current_timestamps = subset[\"time_id\"].to_numpy()\n        \n        # 새로운 타임스탬프 생성\n        new_time_id = np.linspace(current_timestamps.min(), current_timestamps.max(), target_timestamps)\n        \n        # 보간 및 외삽\n        interpolated_values = np.interp(\n            new_time_id,\n            current_timestamps,\n            subset[\"weighted_mean_responder_6\"].to_numpy()\n        )\n        \n        # 결과 저장\n        interpolated_data.append(\n            pd.DataFrame({\n                \"date_id\": [date] * len(new_time_id),\n                \"time_id\": new_time_id,\n                \"weighted_mean_responder_6\": interpolated_values\n            })\n        )\n    return pd.concat(interpolated_data)\n\n# Step 4: 타임스탬프 맞추기\nearly_interpolated = interpolate_and_extrapolate(early_data, 968)\nlate_interpolated = interpolate_and_extrapolate(late_data, 968)\n\n# 데이터 병합 및 정렬\nfinal_interpolated = pd.concat([early_interpolated, late_interpolated])\nfinal_interpolated = final_interpolated.sort_values([\"date_id\", \"time_id\"]).reset_index(drop=True)\n\n# Step 5: 시각화\ndate_choice_early = 500  # date_id <= 676\ndate_choice_late = 700  # date_id > 676\n\nsubset_early = final_interpolated[final_interpolated[\"date_id\"] == date_choice_early]\nsubset_late = final_interpolated[final_interpolated[\"date_id\"] == date_choice_late]\n\n# 그래프 생성\nfig, axes = plt.subplots(1, 2, figsize=(14, 6), sharey=True)\n\n# Early date 시각화\naxes[0].plot(\n    subset_early[\"time_id\"], subset_early[\"weighted_mean_responder_6\"],\n    label=\"Interpolated Data\", linestyle=\"-\", alpha=0.8\n)\naxes[0].set_title(f\"Responder_6 (Date ID {date_choice_early})\")\naxes[0].grid(True, linestyle=\"--\")\n\n# Late date 시각화\naxes[1].plot(\n    subset_late[\"time_id\"], subset_late[\"weighted_mean_responder_6\"],\n    label=\"Interpolated Data\", linestyle=\"-\", alpha=0.8\n)\naxes[1].set_title(f\"Responder_6 (Date ID {date_choice_late})\")\naxes[1].grid(True, linestyle=\"--\")\n\nplt.suptitle(\"Weighted Mean Responder_6 with Interpolation and Extrapolation\")\nplt.show()\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T06:24:20.646092Z","iopub.execute_input":"2024-12-26T06:24:20.646582Z","iopub.status.idle":"2024-12-26T06:24:41.753829Z","shell.execute_reply.started":"2024-12-26T06:24:20.646542Z","shell.execute_reply":"2024-12-26T06:24:41.752358Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"final_interpolated.isna().sum()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T06:24:49.710787Z","iopub.execute_input":"2024-12-26T06:24:49.712097Z","iopub.status.idle":"2024-12-26T06:24:49.733437Z","shell.execute_reply.started":"2024-12-26T06:24:49.712030Z","shell.execute_reply":"2024-12-26T06:24:49.732087Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"final_interpolated.shape","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T06:24:52.661372Z","iopub.execute_input":"2024-12-26T06:24:52.661980Z","iopub.status.idle":"2024-12-26T06:24:52.671426Z","shell.execute_reply.started":"2024-12-26T06:24:52.661923Z","shell.execute_reply":"2024-12-26T06:24:52.669801Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"final_interpolated.head()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T06:24:55.490667Z","iopub.execute_input":"2024-12-26T06:24:55.491112Z","iopub.status.idle":"2024-12-26T06:24:55.503567Z","shell.execute_reply.started":"2024-12-26T06:24:55.491073Z","shell.execute_reply":"2024-12-26T06:24:55.502200Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"final_interpolated.info","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T05:48:49.410053Z","iopub.execute_input":"2024-12-26T05:48:49.410468Z","iopub.status.idle":"2024-12-26T05:48:49.422308Z","shell.execute_reply.started":"2024-12-26T05:48:49.410431Z","shell.execute_reply":"2024-12-26T05:48:49.421217Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport statsmodels.api as sm\nfrom statsmodels.tsa.stattools import adfuller\nimport matplotlib.pyplot as plt\n\n# 1. 메모리 절약형 정상성 분석 함수\ndef test_stationarity_efficient(ts, window=300):\n    # 시계열 샘플링 (메모리 부담 줄이기)\n    sampled_ts = ts[::10]  # 10개 중 1개 샘플링\n\n    # 시각적 검토\n    plt.figure(figsize=(12, 6))\n    plt.plot(sampled_ts, label=\"Sampled Data\", alpha=0.7)\n    plt.title(\"Sampled Responder_6 Time Series\")\n    plt.xlabel(\"Index\")\n    plt.ylabel(\"Responder_6\")\n    plt.legend()\n    plt.grid()\n    plt.show()\n\n    # 이동 평균 및 표준편차 계산\n    rolling_mean = pd.Series(sampled_ts).rolling(window=window, min_periods=1).mean()\n    rolling_std = pd.Series(sampled_ts).rolling(window=window, min_periods=1).std()\n\n    plt.figure(figsize=(12, 6))\n    plt.plot(sampled_ts, label=\"Sampled Data\", alpha=0.7)\n    plt.plot(rolling_mean, label=\"Rolling Mean\", color='orange', alpha=0.9)\n    plt.plot(rolling_std, label=\"Rolling Std\", color='green', alpha=0.9)\n    plt.title(\"Rolling Statistics (Sampled Data)\")\n    plt.legend()\n    plt.grid()\n    plt.show()\n\n    # ADF Test\n    print(\"Results of Augmented Dickey-Fuller Test (Sampled Data):\")\n    result = adfuller(sampled_ts, autolag='AIC')\n    print(f\"ADF Statistic: {result[0]}\")\n    print(f\"p-value: {result[1]}\")\n    print(\"Critical Values:\")\n    for key, value in result[4].items():\n        print(f\"   {key}: {value}\")\n\n    if result[1] <= 0.05:\n        print(\"\\nThe series is stationary (p-value <= 0.05).\")\n    else:\n        print(\"\\nThe series is NOT stationary (p-value > 0.05).\")\n\n# 2. 정상성 분석 실행\n# final_interpolated는 보간 및 외삽이 완료된 pandas DataFrame\nfinal_time_series = final_interpolated[\"weighted_mean_responder_6\"].to_numpy()\n\n# 3. 테스트 실행\ntest_stationarity_efficient(final_time_series)\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T05:48:49.423778Z","iopub.execute_input":"2024-12-26T05:48:49.424159Z","iopub.status.idle":"2024-12-26T05:50:10.968307Z","shell.execute_reply.started":"2024-12-26T05:48:49.424124Z","shell.execute_reply":"2024-12-26T05:50:10.965454Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"final_interpolated.head()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T06:23:30.932628Z","iopub.execute_input":"2024-12-26T06:23:30.933049Z","iopub.status.idle":"2024-12-26T06:23:30.946883Z","shell.execute_reply.started":"2024-12-26T06:23:30.933015Z","shell.execute_reply":"2024-12-26T06:23:30.945368Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"from statsmodels.graphics.tsaplots import plot_acf, plot_pacf\n\n# 1. ACF 플롯\nplt.figure(figsize=(12, 6))\nplot_acf(final_time_series, lags=10, title=\"Autocorrelation Function (ACF)\")\nplt.show()\n\n# 2. PACF 플롯\nplt.figure(figsize=(12, 6))\nplot_pacf(final_time_series, lags=10, title=\"Partial Autocorrelation Function (PACF)\")\nplt.show()\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T05:50:10.970677Z","iopub.execute_input":"2024-12-26T05:50:10.971604Z","iopub.status.idle":"2024-12-26T06:00:36.383313Z","shell.execute_reply.started":"2024-12-26T05:50:10.971542Z","shell.execute_reply":"2024-12-26T06:00:36.381485Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"final_interpolated.head()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T06:23:21.847512Z","iopub.execute_input":"2024-12-26T06:23:21.847977Z","iopub.status.idle":"2024-12-26T06:23:21.860845Z","shell.execute_reply.started":"2024-12-26T06:23:21.847939Z","shell.execute_reply":"2024-12-26T06:23:21.859450Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nfrom statsmodels.tsa.arima.model import ARIMA\nfrom sklearn.metrics import mean_squared_error\n\n# 데이터 로드 및 준비 (final_interpolated 사용)\ndata = final_interpolated['weighted_mean_responder_6'].values\n\n# 학습/검증 데이터 분리\ntrain_size = int(len(data) * 0.8)\ntrain, test = data[:train_size], data[train_size:]\n\n# ARIMA 모델 설정 및 학습\np, d, q = 1, 0, 1  # ARIMA(1, 0, 1)\nmodel = ARIMA(train, order=(p, d, q))\nmodel_fit = model.fit()\n\n# 요약 정보 출력\nprint(model_fit.summary())\n\n# 예측 수행\nforecast_steps = len(test)\nforecast = model_fit.forecast(steps=forecast_steps)\n\n# 시각화\nplt.figure(figsize=(12, 6))\nplt.plot(range(len(train)), train, label='Train Data')\nplt.plot(range(len(train), len(train) + len(test)), test, label='Test Data', color='orange')\nplt.plot(range(len(train), len(train) + len(test)), forecast, label='Forecast', color='green')\nplt.title('ARIMA(1, 0, 1) Forecast')\nplt.legend()\nplt.grid()\nplt.show()\n\n# 성능 평가\nmse = mean_squared_error(test, forecast)\nprint(f'Mean Squared Error (MSE): {mse}')\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T06:00:36.385691Z","iopub.execute_input":"2024-12-26T06:00:36.386191Z","iopub.status.idle":"2024-12-26T06:12:48.967852Z","shell.execute_reply.started":"2024-12-26T06:00:36.386148Z","shell.execute_reply":"2024-12-26T06:12:48.966354Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"final_interpolated.head()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T06:25:12.454557Z","iopub.execute_input":"2024-12-26T06:25:12.455364Z","iopub.status.idle":"2024-12-26T06:25:12.469400Z","shell.execute_reply.started":"2024-12-26T06:25:12.455306Z","shell.execute_reply":"2024-12-26T06:25:12.467889Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nfrom statsmodels.tsa.statespace.sarimax import SARIMAX\nfrom sklearn.metrics import mean_squared_error\nimport matplotlib.pyplot as plt\n\n# 데이터 로드 (예제)\n# final_interpolated는 이미 보간된 데이터\nfinal_interpolated = final_interpolated.sort_values(by=\"time_id\")  # time_id 기준 정렬\n\n# time_id를 인덱스로 설정\nfinal_interpolated.set_index(\"time_id\", inplace=True)\nprint(final_interpolated.head())\n# 훈련 데이터와 테스트 데이터 분리\ntrain_data = final_interpolated.iloc[:int(0.8 * len(final_interpolated))]\ntest_data = final_interpolated.iloc[int(0.8 * len(final_interpolated)):]\n\n# 응답 변수 (responder_6)\ntrain_y = train_data[\"weighted_mean_responder_6\"]\ntest_y = test_data[\"weighted_mean_responder_6\"]\n\n# SARIMAX 모델 학습\norder = (1, 0, 1)  # ARIMA(p, d, q)\nseasonal_order = (0, 0, 0, 0)  # SARIMA(P, D, Q, S)\n\nmodel = SARIMAX(\n    train_y,\n    order=order,\n    seasonal_order=seasonal_order,\n    enforce_stationarity=False,\n    enforce_invertibility=False,\n)\nmodel_fit = model.fit(disp=False)\n\n# 테스트 데이터 예측\nforecast_steps = len(test_y)\nforecast = model_fit.forecast(steps=forecast_steps)\n\n# MSE 평가\nmse = mean_squared_error(test_y, forecast)\nprint(f\"Mean Squared Error (MSE): {mse}\")\n\n# 시각화\nplt.figure(figsize=(12, 6))\nplt.plot(train_data.index, train_y, label=\"Train Data\")\nplt.plot(test_data.index, test_y, label=\"Test Data\", color=\"orange\")\nplt.plot(test_data.index, forecast, label=\"Forecast\", color=\"green\")\nplt.title(f\"SARIMAX({order}) Forecast\")\nplt.legend()\nplt.grid()\nplt.show()\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-12-26T06:25:22.440442Z","iopub.execute_input":"2024-12-26T06:25:22.440987Z","iopub.status.idle":"2024-12-26T06:29:38.384838Z","shell.execute_reply.started":"2024-12-26T06:25:22.440940Z","shell.execute_reply":"2024-12-26T06:29:38.382790Z"}},"outputs":[],"execution_count":null}]}