{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.12.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"gpu","dataSources":[{"sourceId":11000,"databundleVersionId":875412,"sourceType":"competition"}],"dockerImageVersionId":31236,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"# ============================================================\n# LANL Earthquake Prediction\n# Final Lightweight LSTM-only Ensemble (Private-oriented)\n# ============================================================\n\n# ===============================\n# 0. Imports & Device\n# ===============================\nimport os, gc, random\nimport numpy as np\nimport pandas as pd\nimport torch\nimport torch.nn as nn\nimport pywt\nfrom sklearn.preprocessing import StandardScaler\nfrom scipy.stats import skew\nfrom torch.utils.data import Dataset, DataLoader\n\ndevice = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\nprint(\"device:\", device)\n\n# ===============================\n# 1. Reproducibility\n# ===============================\ndef seed_everything(seed):\n    random.seed(seed)\n    np.random.seed(seed)\n    torch.manual_seed(seed)\n    torch.cuda.manual_seed_all(seed)\n\n# ===============================\n# 2. Signal Preprocess (超重要)\n# ===============================\ndef preprocess_signal(x):\n    x = x + np.random.normal(0, 0.5, len(x))\n    x = x - np.median(x)\n    return x.astype(np.float32)\n\n# ===============================\n# 3. Patch Generator\n# ===============================\ndef patch_generator(csv_path, patch_size, overlap, chunksize=1_000_000):\n    step = patch_size - overlap\n    buf_x = np.empty(0, dtype=np.float32)\n    buf_y = np.empty(0, dtype=np.float32)\n\n    for chunk in pd.read_csv(csv_path, chunksize=chunksize):\n        x = preprocess_signal(chunk.acoustic_data.values)\n        y = chunk.time_to_failure.values.astype(np.float32)\n\n        buf_x = np.concatenate([buf_x, x])\n        buf_y = np.concatenate([buf_y, y])\n\n        while len(buf_x) >= patch_size:\n            yield buf_x[:patch_size], buf_y[patch_size - 1]\n            buf_x = buf_x[step:]\n            buf_y = buf_y[step:]\n\n        del chunk, x, y\n        gc.collect()\n\n# ===============================\n# 4. Feature Engineering（最軽量）\n# ===============================\ndef fft_low(x):\n    p = np.abs(np.fft.rfft(x))\n    return np.log1p(p[:20].mean())\n\ndef wavelet_energy(x):\n    _, cD = pywt.dwt(x, \"db1\")\n    return np.log1p(np.mean(np.abs(cD)))\n\ndef block_features(x):\n    q = np.percentile(x, [25, 50, 75])\n    return np.array([\n        x.mean(),\n        x.std(),\n        skew(x),\n        np.mean(np.abs(x)),\n        np.sqrt(np.mean(x ** 2)),\n        q[0], q[1], q[2],\n        q[2] - q[0],\n        fft_low(x),\n        wavelet_energy(x)\n    ], dtype=np.float32)\n\nFEATURE_DIM = 11\n\ndef sequence_features(x, n_blocks=80):\n    blocks = np.array_split(x, n_blocks)\n    return np.stack([block_features(b) for b in blocks])\n\n# ===============================\n# 5. Dataset\n# ===============================\nclass LANLDataset(Dataset):\n    def __init__(self, csv_path, patch_size, overlap, limit):\n        self.samples = []\n        for i, (x, y) in enumerate(patch_generator(csv_path, patch_size, overlap)):\n            self.samples.append((sequence_features(x), y))\n            if i >= limit:\n                break\n\n    def __len__(self):\n        return len(self.samples)\n\n    def __getitem__(self, idx):\n        x, y = self.samples[idx]\n        return torch.tensor(x), torch.tensor(y)\n\n# ===============================\n# 6. Weak Attention Pooling\n# ===============================\nclass WeakAttentionPooling(nn.Module):\n    def __init__(self, h):\n        super().__init__()\n        self.attn = nn.Linear(h, 1)\n\n    def forward(self, x):\n        w = torch.tanh(self.attn(x)).squeeze(-1)\n        w = w / (w.abs().mean(dim=1, keepdim=True) + 1e-6)\n        return torch.mean(x * w.unsqueeze(-1), dim=1)\n\n# ===============================\n# 7. LSTM Model（唯一のモデル）\n# ===============================\nclass LSTMModel(nn.Module):\n    def __init__(self):\n        super().__init__()\n        self.lstm = nn.LSTM(\n            FEATURE_DIM,\n            128,\n            num_layers=2,\n            batch_first=True,\n            bidirectional=True,\n            dropout=0.5\n        )\n        self.pool = WeakAttentionPooling(256)\n        self.fc = nn.Linear(256, 1)\n\n    def forward(self, x):\n        h, _ = self.lstm(x)\n        return self.fc(self.pool(h)).squeeze()\n\n# ===============================\n# 8. Train + Predict\n# ===============================\ndef train_and_predict(root, overlap_ratio, seed):\n    seed_everything(seed)\n\n    patch_size = 150_000\n    overlap = int(patch_size * overlap_ratio)\n    limit = 1800   # Private寄り\n\n    dataset = LANLDataset(root + \"train.csv\", patch_size, overlap, limit)\n\n    X_all = np.vstack([s[0].reshape(-1, FEATURE_DIM) for s in dataset.samples])\n    scaler = StandardScaler().fit(X_all)\n\n    for i in range(len(dataset.samples)):\n        x, y = dataset.samples[i]\n        dataset.samples[i] = (scaler.transform(x), y)\n\n    split = int(len(dataset) * 0.85)\n    train_ds = torch.utils.data.Subset(dataset, range(split))\n    valid_ds = torch.utils.data.Subset(dataset, range(split, len(dataset)))\n\n    train_loader = DataLoader(train_ds, batch_size=32, shuffle=True)\n    valid_loader = DataLoader(valid_ds, batch_size=64)\n\n    model = LSTMModel().to(device)\n    opt = torch.optim.AdamW(model.parameters(), lr=3e-4)\n    sched = torch.optim.lr_scheduler.ReduceLROnPlateau(opt, patience=3)\n    loss_fn = nn.SmoothL1Loss(beta=0.5)\n\n    for epoch in range(20):\n        model.train()\n        for xb, yb in train_loader:\n            xb, yb = xb.to(device), yb.to(device)\n            opt.zero_grad()\n            loss_fn(model(xb), yb).backward()\n            opt.step()\n\n        model.eval()\n        va = np.mean([\n            loss_fn(model(xb.to(device)), yb.to(device)).item()\n            for xb, yb in valid_loader\n        ])\n        sched.step(va)\n\n    # ---- Test ----\n    X_test = []\n    for f in sorted(os.listdir(root + \"test\")):\n        x = preprocess_signal(\n            pd.read_csv(root + \"test/\" + f).acoustic_data.values\n        )\n        X_test.append(scaler.transform(sequence_features(x)))\n\n    X_test = torch.tensor(np.array(X_test)).to(device)\n    with torch.no_grad():\n        preds = model(X_test).cpu().numpy()\n\n    return preds\n\n# ===============================\n# 9. Ensemble（Public犠牲）\n# ===============================\nroot = \"../input/LANL-Earthquake-Prediction/\"\noverlaps = [0.15,0.20]\nseeds = [42,128,343,777]\n\nall_preds = []\nfor ov in overlaps:\n    for sd in seeds:\n        print(f\"Running overlap={ov}, seed={sd}\")\n        all_preds.append(train_and_predict(root, ov, sd))\n\nfinal_pred = np.median(np.stack(all_preds), axis=0)\n\n# ===============================\n# 10. Submission\n# ===============================\nsample = pd.read_csv(root + \"sample_submission.csv\")\nsample[\"time_to_failure\"] = final_pred\nsample.to_csv(\"submission_final_private_push.csv\", index=False)\n\nprint(\"DONE\")\n","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2026-01-15T02:24:18.039906Z","iopub.execute_input":"2026-01-15T02:24:18.040558Z","iopub.status.idle":"2026-01-15T03:24:14.051579Z","shell.execute_reply.started":"2026-01-15T02:24:18.040525Z","shell.execute_reply":"2026-01-15T03:24:14.050763Z"}},"outputs":[],"execution_count":null}]}