{"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":"nvidiaTeslaT4","dataSources":[{"sourceType":"competition","sourceId":35297,"databundleVersionId":3589144}],"dockerImageVersionId":31260,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python Docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load\n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\n\n# Input data files are available in the read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\nimport os\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    for filename in filenames:\n        print(os.path.join(dirname, filename))\n\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2026-03-03T23:06:11.490241Z","iopub.execute_input":"2026-03-03T23:06:11.491049Z","iopub.status.idle":"2026-03-03T23:06:11.496971Z","shell.execute_reply.started":"2026-03-03T23:06:11.491013Z","shell.execute_reply":"2026-03-03T23:06:11.496158Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ==============================\n# Kaggle Environment Imports\n# ==============================\nimport os\nimport warnings\nwarnings.filterwarnings('ignore')\n\nimport numpy as np\nimport pandas as pd\n\nimport torch\nfrom torch import nn, optim\nfrom torch.utils.data import DataLoader, TensorDataset\n\nfrom sklearn.preprocessing import LabelEncoder, RobustScaler\nfrom sklearn.metrics import roc_auc_score, average_precision_score\nfrom sklearn.model_selection import KFold\n\nfrom scipy.signal import savgol_filter\n\n# ==============================\n# 1. PHYSICS-INFORMED DATA ENGINE\n# ==============================\ndef prepare_physics_sota_data(path):\n    print(\"🧬 Engineering Thermal Inertia & Kinetic Momentum Manifolds...\")\n\n    df = pd.read_csv(path)\n\n    # Kaggle-safe timestamp parsing\n    df['timestamp'] = pd.to_datetime(df['timestamp'], errors='coerce')\n    df = df.dropna(subset=['timestamp'])\n\n    df = df.sort_values(['building_id', 'timestamp'])\n\n    # --- Thermodynamic baseline ---\n    df['meter_log'] = np.log1p(df['meter_reading'])\n    df['hour'] = df['timestamp'].dt.hour\n    df['dow'] = df['timestamp'].dt.dayofweek\n\n    norm_stats = (\n        df[df['anomaly'] == 0]\n        .groupby(['building_id', 'dow', 'hour'])['meter_log']\n        .agg(['median', 'std'])\n        .reset_index()\n    )\n    norm_stats.columns = ['building_id', 'dow', 'hour', 'h_med', 'h_std']\n    df = df.merge(norm_stats, on=['building_id', 'dow', 'hour'], how='left')\n\n    df['h_med'] = df['h_med'].fillna(\n        df.groupby('building_id')['meter_log'].transform('median')\n    )\n    df['h_std'] = df['h_std'].fillna(\n        df.groupby('building_id')['meter_log'].transform('std')\n    ).replace(0, 1.0)\n\n    # --- Physics Feature 1: Thermal Inertia ---\n    df['temp_inertia'] = df.groupby('building_id')['air_temperature'] \\\n        .transform(lambda x: x.rolling(6, min_periods=1).mean())\n    df['thermal_p'] = df['air_temperature'] - df['temp_inertia']\n\n    # --- Physics Feature 2: Kinetic Shock ---\n    df['shock'] = df['meter_log'] - df['h_med']\n    df['shock_v'] = df.groupby('building_id')['shock'].diff().fillna(0)\n    df['kinetic_s'] = df.groupby('building_id')['shock_v'].diff().fillna(0)\n\n    # --- Robust Scaling ---\n    scaler = RobustScaler()\n    df['shock_scaled'] = df.groupby('building_id')['shock'].transform(\n        lambda x: scaler.fit_transform(x.values.reshape(-1, 1)).flatten()\n        if x.std() > 0 else x\n    )\n\n    # --- Encoding ---\n    le_b, le_s, le_u = LabelEncoder(), LabelEncoder(), LabelEncoder()\n    df['b_idx'] = le_b.fit_transform(df['building_id'])\n    df['s_idx'] = le_s.fit_transform(df['site_id'])\n    df['u_idx'] = le_u.fit_transform(df['primary_use'])\n\n    feat_cols = ['shock_scaled', 'thermal_p', 'kinetic_s', 'air_temperature']\n    df = df.replace([np.inf, -np.inf], np.nan).fillna(0)\n\n    # --- Same research split ---\n    train_df = df[df['anomaly'] == 0].sample(frac=0.45, random_state=42)\n    test_df = df.drop(train_df.index).sort_values(['building_id', 'timestamp'])\n\n    return {\n        'train': (\n            torch.FloatTensor(train_df[feat_cols].values),\n            torch.LongTensor(train_df['b_idx'].values),\n            torch.LongTensor(train_df['s_idx'].values),\n            torch.LongTensor(train_df['u_idx'].values)\n        ),\n        'test': (\n            torch.FloatTensor(test_df[feat_cols].values),\n            torch.LongTensor(test_df['b_idx'].values),\n            torch.LongTensor(test_df['s_idx'].values),\n            torch.LongTensor(test_df['u_idx'].values),\n            test_df['anomaly'].values,\n            test_df['building_id'].values\n        ),\n        'meta': {\n            'num_b': len(le_b.classes_),\n            'num_s': len(le_s.classes_),\n            'num_u': len(le_u.classes_),\n            'feats': feat_cols\n        }\n    }\n\n# ==============================\n# 2. MODEL ARCHITECTURE (PIGD)\n# ==============================\nclass PIGD_Block(nn.Module):\n    def __init__(self, dim, ctx_dim):\n        super().__init__()\n        self.ln = nn.LayerNorm(dim)\n        self.film_gamma = nn.Linear(ctx_dim, dim)\n        self.film_beta = nn.Linear(ctx_dim, dim)\n        self.gate = nn.Sequential(nn.Linear(dim, dim), nn.Sigmoid())\n        self.net = nn.Sequential(\n            nn.Linear(dim, dim * 2),\n            nn.SiLU(),\n            nn.Linear(dim * 2, dim)\n        )\n\n    def forward(self, x, ctx):\n        h = self.ln(x)\n        h = h * (1 + self.film_gamma(ctx)) + self.film_beta(ctx)\n        return x + self.gate(h) * self.net(h)\n\nclass PIGD_Model(nn.Module):\n    def __init__(self, input_dim, num_b, num_s, num_u):\n        super().__init__()\n        self.b_emb = nn.Embedding(num_b, 64)\n        self.s_emb = nn.Embedding(num_s, 16)\n        self.u_emb = nn.Embedding(num_u, 16)\n        self.t_emb = nn.Sequential(nn.Linear(1, 32), nn.SiLU())\n\n        ctx_dim = 64 + 16 + 16 + 32\n        self.proj_in = nn.Linear(input_dim, 128)\n        self.layers = nn.ModuleList([PIGD_Block(128, ctx_dim) for _ in range(4)])\n        self.proj_out = nn.Linear(128, input_dim)\n\n    def forward(self, x, t, b, s, u):\n        te = self.t_emb(t.view(-1, 1).float() / 100.0)\n        ctx = torch.cat([self.b_emb(b), self.s_emb(s), self.u_emb(u), te], dim=-1)\n        h = self.proj_in(x)\n        for layer in self.layers:\n            h = layer(h, ctx)\n        return self.proj_out(h)\n\n# ==============================\n# 3. EXECUTION WITH K-FOLD = 3\n# ==============================\ndef run_piml_engine(data):\n    device = 'cuda' if torch.cuda.is_available() else 'cpu'\n    meta = data['meta']\n\n    x_all, b_all, s_all, u_all = data['train']\n    kf = KFold(n_splits=3, shuffle=True, random_state=42)\n\n    aucs, prs = [], []\n\n    for fold, (idx, _) in enumerate(kf.split(x_all)):\n        print(f\"\\n🔥 Fold {fold+1}/3 Training...\")\n\n        model = PIGD_Model(\n            len(meta['feats']), meta['num_b'], meta['num_s'], meta['num_u']\n        ).to(device)\n\n        optimizer = optim.AdamW(model.parameters(), lr=1e-3, weight_decay=1e-5)\n        alphas = torch.cumprod(1. - torch.linspace(1e-4, 0.02, 100), 0).to(device)\n\n        loader = DataLoader(\n            TensorDataset(x_all[idx], b_all[idx], s_all[idx], u_all[idx]),\n            batch_size=1024,\n            shuffle=True\n        )\n\n        for epoch in range(15):\n            model.train()\n            for x0, b, s, u in loader:\n                x0, b, s, u = x0.to(device), b.to(device), s.to(device), u.to(device)\n                t = torch.randint(0, 100, (x0.size(0),), device=device)\n                noise = torch.randn_like(x0)\n                a_t = alphas[t].view(-1, 1)\n                xt = torch.sqrt(a_t)*x0 + torch.sqrt(1-a_t)*noise\n                loss = nn.MSELoss()(model(xt, t, b, s, u), x0)\n                optimizer.zero_grad()\n                loss.backward()\n                optimizer.step()\n\n        model.eval()\n        x_te, b_te, s_te, u_te, y_te, b_ids = data['test']\n        x_te, b_te, s_te, u_te = x_te.to(device), b_te.to(device), s_te.to(device), u_te.to(device)\n\n        with torch.no_grad():\n            t = torch.ones(x_te.size(0), device=device).long() * 25\n            a_t = alphas[t].view(-1, 1)\n            xt = torch.sqrt(a_t)*x_te + torch.sqrt(1-a_t)*torch.randn_like(x_te)\n            recon = model(xt, t, b_te, s_te, u_te)\n\n            err = torch.abs(x_te[:, 0] - recon[:, 0]) + 0.5 * torch.abs(x_te[:, 2] - recon[:, 2])\n            raw_scores = err.cpu().numpy()\n\n            df_err = pd.DataFrame({'err': raw_scores, 'b': b_ids})\n            scores = df_err.groupby('b')['err'].transform(\n                lambda x: savgol_filter(x, min(len(x), 11), 2) if len(x) > 11 else x\n            )\n\n        auc = roc_auc_score(y_te, scores)\n        pr = average_precision_score(y_te, scores)\n\n        print(f\"Fold {fold+1} ROC-AUC: {auc:.4f}\")\n        print(f\"Fold {fold+1} PR-AUC : {pr:.4f}\")\n\n        aucs.append(auc)\n        prs.append(pr)\n\n    print(\"\\n🏆 FINAL K-FOLD (k=3) RESULTS\")\n    print(f\"Mean ROC-AUC: {np.mean(aucs):.4f}\")\n    print(f\"Mean PR-AUC : {np.mean(prs):.4f}\")\n\n# ==============================\n# MAIN (KAGGLE PATH)\n# ==============================\nif __name__ == \"__main__\":\n    path = \"/kaggle/input/energy-anomaly-detection/train_features.csv\"\n    if os.path.exists(path):\n        data = prepare_physics_sota_data(path)\n        run_piml_engine(data)\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-28T17:52:41.016912Z","iopub.execute_input":"2026-01-28T17:52:41.017203Z","iopub.status.idle":"2026-01-28T18:00:42.007904Z","shell.execute_reply.started":"2026-01-28T17:52:41.017184Z","shell.execute_reply":"2026-01-28T18:00:42.007153Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# =========================================================\n# Kaggle Environment & Imports\n# =========================================================\nimport os\nimport warnings\nwarnings.filterwarnings(\"ignore\")\n\nimport numpy as np\nimport pandas as pd\n\nimport torch\nfrom torch import nn, optim\nfrom torch.utils.data import DataLoader, TensorDataset\n\nfrom sklearn.preprocessing import LabelEncoder\nfrom sklearn.metrics import roc_auc_score, average_precision_score\nfrom sklearn.model_selection import KFold\n\nfrom scipy.signal import savgol_filter\n\n# =========================================================\n# DEVICE\n# =========================================================\ndevice = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\nprint(\"Device:\", device)\nif device.type == \"cuda\":\n    print(\"GPU:\", torch.cuda.get_device_name(0))\n\n# =========================================================\n# 1. DATA PREPARATION (Original logic preserved)\n# =========================================================\ndef prepare_physics_sota_data(path):\n    print(\"🧬 Engineering Thermal Inertia & Kinetic Momentum Manifolds...\")\n\n    df = pd.read_csv(path)\n\n    # ✅ FIXED TIMESTAMP PARSING (Kaggle-safe)\n    df[\"timestamp\"] = pd.to_datetime(df[\"timestamp\"], format=\"mixed\", errors=\"coerce\")\n    df = df.dropna(subset=[\"timestamp\"])\n\n    df = df.sort_values([\"building_id\", \"timestamp\"])\n\n    # --- Thermodynamic baseline ---\n    df[\"meter_log\"] = np.log1p(df[\"meter_reading\"])\n    df[\"hour\"] = df[\"timestamp\"].dt.hour\n    df[\"dow\"] = df[\"timestamp\"].dt.dayofweek\n\n    norm = (\n        df[df[\"anomaly\"] == 0]\n        .groupby([\"building_id\", \"dow\", \"hour\"])[\"meter_log\"]\n        .agg([\"median\", \"std\"])\n        .reset_index()\n    )\n    norm.columns = [\"building_id\", \"dow\", \"hour\", \"h_med\", \"h_std\"]\n    df = df.merge(norm, on=[\"building_id\", \"dow\", \"hour\"], how=\"left\")\n\n    df[\"h_med\"] = df[\"h_med\"].fillna(\n        df.groupby(\"building_id\")[\"meter_log\"].transform(\"median\")\n    )\n    df[\"h_std\"] = df[\"h_std\"].replace(0, 1).fillna(1)\n\n    # --- Physics Feature 1: Thermal Inertia ---\n    df[\"temp_inertia\"] = (\n        df.groupby(\"building_id\")[\"air_temperature\"]\n        .transform(lambda x: x.rolling(6, min_periods=1).mean())\n    )\n    df[\"thermal_p\"] = df[\"air_temperature\"] - df[\"temp_inertia\"]\n\n    # --- Physics Feature 2: Kinetic Shock ---\n    df[\"shock\"] = df[\"meter_log\"] - df[\"h_med\"]\n    df[\"shock_v\"] = df.groupby(\"building_id\")[\"shock\"].diff().fillna(0)\n    df[\"kinetic_s\"] = df.groupby(\"building_id\")[\"shock_v\"].diff().fillna(0)\n\n    # --- Robust scaling (same idea, safe) ---\n    df[\"shock_scaled\"] = df.groupby(\"building_id\")[\"shock\"].transform(\n        lambda x: (x - x.median()) /\n        (x.quantile(0.75) - x.quantile(0.25) + 1e-6)\n    )\n\n    feat_cols = [\"shock_scaled\", \"thermal_p\", \"kinetic_s\", \"air_temperature\"]\n\n    # Stability (NO aggressive clipping)\n    df[feat_cols] = df[feat_cols].replace([np.inf, -np.inf], 0).fillna(0)\n    df[feat_cols] = df[feat_cols].clip(-25, 25)\n\n    # --- Encoding ---\n    le_b, le_s, le_u = LabelEncoder(), LabelEncoder(), LabelEncoder()\n    df[\"b_idx\"] = le_b.fit_transform(df[\"building_id\"])\n    df[\"s_idx\"] = le_s.fit_transform(df[\"site_id\"])\n    df[\"u_idx\"] = le_u.fit_transform(df[\"primary_use\"])\n\n    # --- Same research split ---\n    train_df = df[df[\"anomaly\"] == 0].sample(frac=0.45, random_state=42)\n    test_df = df.drop(train_df.index).sort_values([\"building_id\", \"timestamp\"])\n\n    return {\n        \"train\": (\n            torch.tensor(train_df[feat_cols].values, dtype=torch.float32),\n            torch.tensor(train_df[\"b_idx\"].values),\n            torch.tensor(train_df[\"s_idx\"].values),\n            torch.tensor(train_df[\"u_idx\"].values),\n        ),\n        \"test\": (\n            torch.tensor(test_df[feat_cols].values, dtype=torch.float32),\n            torch.tensor(test_df[\"b_idx\"].values),\n            torch.tensor(test_df[\"s_idx\"].values),\n            torch.tensor(test_df[\"u_idx\"].values),\n            test_df[\"anomaly\"].values.astype(np.int32),\n            test_df[\"building_id\"].values,\n        ),\n        \"meta\": {\n            \"num_b\": len(le_b.classes_),\n            \"num_s\": len(le_s.classes_),\n            \"num_u\": len(le_u.classes_),\n            \"feats\": feat_cols,\n        },\n    }\n\n# =========================================================\n# 2. MODEL (UNCHANGED CORE LOGIC)\n# =========================================================\nclass PIGD_Block(nn.Module):\n    def __init__(self, d, c):\n        super().__init__()\n        self.ln = nn.LayerNorm(d)\n        self.g = nn.Linear(c, d)\n        self.b = nn.Linear(c, d)\n        self.gate = nn.Sequential(nn.Linear(d, d), nn.Sigmoid())\n        self.net = nn.Sequential(\n            nn.Linear(d, d * 2),\n            nn.SiLU(),\n            nn.Linear(d * 2, d)\n        )\n\n    def forward(self, x, ctx):\n        h = self.ln(x)\n        h = h * (1 + self.g(ctx)) + self.b(ctx)\n        return x + self.gate(h) * self.net(h)\n\nclass PIGD_Model(nn.Module):\n    def __init__(self, inp, nb, ns, nu):\n        super().__init__()\n        self.b = nn.Embedding(nb, 64)\n        self.s = nn.Embedding(ns, 16)\n        self.u = nn.Embedding(nu, 16)\n        self.t = nn.Sequential(nn.Linear(1, 32), nn.SiLU())\n\n        ctx_dim = 64 + 16 + 16 + 32\n        self.inp = nn.Linear(inp, 256)\n        self.blocks = nn.ModuleList([PIGD_Block(256, ctx_dim) for _ in range(6)])\n        self.out = nn.Linear(256, inp)\n\n    def forward(self, x, t, b, s, u):\n        te = self.t(t.float().unsqueeze(1) / 100)\n        ctx = torch.cat([self.b(b), self.s(s), self.u(u), te], dim=1)\n        h = self.inp(x)\n        for blk in self.blocks:\n            h = blk(h, ctx)\n        return self.out(h)\n\n# =========================================================\n# SAFE SAVITZKY–GOLAY\n# =========================================================\ndef safe_savgol(x):\n    x = np.asarray(x, dtype=np.float32)\n    if len(x) < 7 or np.allclose(x, x[0]):\n        return x\n    w = min(11, len(x) // 2 * 2 + 1)\n    try:\n        return savgol_filter(x, w, 2, mode=\"interp\")\n    except:\n        return x\n\n# =========================================================\n# 3. EXECUTION WITH K-FOLD = 3 (NON-MISDIRECTING)\n# =========================================================\ndef run_piml_engine(data):\n    meta = data[\"meta\"]\n    x_all, b_all, s_all, u_all = data[\"train\"]\n\n    kf = KFold(n_splits=3, shuffle=True, random_state=42)\n\n    aucs, prs = [], []\n\n    for fold, (idx, _) in enumerate(kf.split(x_all)):\n        print(f\"\\n🔥 Fold {fold+1}/3 Training...\")\n\n        model = PIGD_Model(\n            len(meta[\"feats\"]),\n            meta[\"num_b\"],\n            meta[\"num_s\"],\n            meta[\"num_u\"]\n        ).to(device)\n\n        opt = optim.AdamW(model.parameters(), lr=1e-3)\n        scaler = torch.cuda.amp.GradScaler()\n        alphas = torch.cumprod(\n            1 - torch.linspace(1e-4, 0.02, 100, device=device), 0\n        )\n\n        loader = DataLoader(\n            TensorDataset(x_all[idx], b_all[idx], s_all[idx], u_all[idx]),\n            batch_size=2048,\n            shuffle=True\n        )\n\n        # ---- Training (same logic, repeated) ----\n        model.train()\n        for e in range(18):\n            for x0, b0, s0, u0 in loader:\n                x0, b0, s0, u0 = (\n                    x0.to(device),\n                    b0.to(device),\n                    s0.to(device),\n                    u0.to(device),\n                )\n\n                t = torch.randint(0, 100, (x0.size(0),), device=device)\n                a = alphas[t].unsqueeze(1)\n                xt = torch.sqrt(a) * x0 + torch.sqrt(1 - a) * torch.randn_like(x0)\n\n                opt.zero_grad()\n                with torch.cuda.amp.autocast():\n                    loss = nn.MSELoss()(model(xt, t, b0, s0, u0), x0)\n\n                scaler.scale(loss).backward()\n                scaler.step(opt)\n                scaler.update()\n\n            print(f\"Epoch {e+1} | Loss {loss.item():.4f}\")\n\n        # ---- Inference (unchanged) ----\n        model.eval()\n        x, b, s, u, y, bid = data[\"test\"]\n        x, b, s, u = x.to(device), b.to(device), s.to(device), u.to(device)\n\n        with torch.no_grad():\n            t = torch.full((x.size(0),), 25, device=device)\n            a = alphas[t].unsqueeze(1)\n            xt = torch.sqrt(a) * x + torch.sqrt(1 - a) * torch.randn_like(x)\n            r = model(xt, t, b, s, u)\n\n            err = (\n                torch.abs(x[:, 0] - r[:, 0]) +\n                0.75 * torch.abs(x[:, 2] - r[:, 2]) +\n                0.5 * torch.abs(x[:, 1] - r[:, 1])\n            ).cpu().numpy()\n\n        df = pd.DataFrame({\"e\": err, \"b\": bid})\n        smooth = df.groupby(\"b\")[\"e\"].transform(safe_savgol)\n        smooth = np.nan_to_num(smooth, 0)\n\n        auc = roc_auc_score(y, smooth)\n        pr = average_precision_score(y, smooth)\n\n        print(f\"Fold {fold+1} ROC-AUC: {auc:.4f}\")\n        print(f\"Fold {fold+1} PR-AUC : {pr:.4f}\")\n\n        aucs.append(auc)\n        prs.append(pr)\n\n    print(\"\\n🏆 FINAL K-FOLD (k=3) RESULTS\")\n    print(f\"Mean ROC-AUC: {np.mean(aucs):.4f}\")\n    print(f\"Mean PR-AUC : {np.mean(prs):.4f}\")\n\n# =========================================================\n# RUN (KAGGLE PATH)\n# =========================================================\nif __name__ == \"__main__\":\n    path = \"/kaggle/input/energy-anomaly-detection/train_features.csv\"\n    if os.path.exists(path):\n        data = prepare_physics_sota_data(path)\n        run_piml_engine(data)\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-28T18:06:57.676421Z","iopub.execute_input":"2026-01-28T18:06:57.67701Z","iopub.status.idle":"2026-01-28T18:20:14.014315Z","shell.execute_reply.started":"2026-01-28T18:06:57.676979Z","shell.execute_reply":"2026-01-28T18:20:14.013218Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# =========================================================\n# Kaggle Environment & Imports\n# =========================================================\nimport os\nimport warnings\nwarnings.filterwarnings(\"ignore\")\n\nimport numpy as np\nimport pandas as pd\n\nimport torch\nfrom torch import nn, optim\nfrom torch.utils.data import DataLoader, TensorDataset\n\nfrom sklearn.preprocessing import LabelEncoder\nfrom sklearn.metrics import roc_auc_score, average_precision_score\nfrom sklearn.model_selection import KFold\n\nfrom scipy.signal import savgol_filter\n\n# =========================================================\n# DEVICE\n# =========================================================\ndevice = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\nprint(\"Device:\", device)\nif device.type == \"cuda\":\n    print(\"GPU:\", torch.cuda.get_device_name(0))\n\n# =========================================================\n# 1. DATA PREPARATION (ORIGINAL LOGIC PRESERVED)\n# =========================================================\ndef prepare_physics_sota_data(path):\n    print(\"🧬 Engineering Thermal Inertia & Kinetic Momentum Manifolds...\")\n\n    df = pd.read_csv(path)\n\n    # Fixed timestamp parsing\n    df[\"timestamp\"] = pd.to_datetime(df[\"timestamp\"], format=\"mixed\", errors=\"coerce\")\n    df = df.dropna(subset=[\"timestamp\"])\n    df = df.sort_values([\"building_id\", \"timestamp\"])\n\n    # Thermodynamic baseline\n    df[\"meter_log\"] = np.log1p(df[\"meter_reading\"])\n    df[\"hour\"] = df[\"timestamp\"].dt.hour\n    df[\"dow\"] = df[\"timestamp\"].dt.dayofweek\n\n    norm = (\n        df[df[\"anomaly\"] == 0]\n        .groupby([\"building_id\", \"dow\", \"hour\"])[\"meter_log\"]\n        .agg([\"median\", \"std\"])\n        .reset_index()\n    )\n    norm.columns = [\"building_id\", \"dow\", \"hour\", \"h_med\", \"h_std\"]\n    df = df.merge(norm, on=[\"building_id\", \"dow\", \"hour\"], how=\"left\")\n\n    df[\"h_med\"] = df[\"h_med\"].fillna(\n        df.groupby(\"building_id\")[\"meter_log\"].transform(\"median\")\n    )\n    df[\"h_std\"] = df[\"h_std\"].replace(0, 1).fillna(1)\n\n    # Physics features\n    df[\"temp_inertia\"] = (\n        df.groupby(\"building_id\")[\"air_temperature\"]\n        .transform(lambda x: x.rolling(6, min_periods=1).mean())\n    )\n    df[\"thermal_p\"] = df[\"air_temperature\"] - df[\"temp_inertia\"]\n\n    df[\"shock\"] = df[\"meter_log\"] - df[\"h_med\"]\n    df[\"shock_v\"] = df.groupby(\"building_id\")[\"shock\"].diff().fillna(0)\n    df[\"kinetic_s\"] = df.groupby(\"building_id\")[\"shock_v\"].diff().fillna(0)\n\n    # Robust scaling\n    df[\"shock_scaled\"] = df.groupby(\"building_id\")[\"shock\"].transform(\n        lambda x: (x - x.median()) /\n        (x.quantile(0.75) - x.quantile(0.25) + 1e-6)\n    )\n\n    feat_cols = [\"shock_scaled\", \"thermal_p\", \"kinetic_s\", \"air_temperature\"]\n\n    df[feat_cols] = df[feat_cols].replace([np.inf, -np.inf], 0).fillna(0)\n    df[feat_cols] = df[feat_cols].clip(-25, 25)\n\n    # Encoding\n    le_b, le_s, le_u = LabelEncoder(), LabelEncoder(), LabelEncoder()\n    df[\"b_idx\"] = le_b.fit_transform(df[\"building_id\"])\n    df[\"s_idx\"] = le_s.fit_transform(df[\"site_id\"])\n    df[\"u_idx\"] = le_u.fit_transform(df[\"primary_use\"])\n\n    train_df = df[df[\"anomaly\"] == 0].sample(frac=0.45, random_state=42)\n    test_df = df.drop(train_df.index).sort_values([\"building_id\", \"timestamp\"])\n\n    return {\n        \"train\": (\n            torch.tensor(train_df[feat_cols].values, dtype=torch.float32),\n            torch.tensor(train_df[\"b_idx\"].values),\n            torch.tensor(train_df[\"s_idx\"].values),\n            torch.tensor(train_df[\"u_idx\"].values),\n        ),\n        \"test\": (\n            torch.tensor(test_df[feat_cols].values, dtype=torch.float32),\n            torch.tensor(test_df[\"b_idx\"].values),\n            torch.tensor(test_df[\"s_idx\"].values),\n            torch.tensor(test_df[\"u_idx\"].values),\n            test_df[\"anomaly\"].values.astype(np.int32),\n            test_df[\"building_id\"].values,\n        ),\n        \"meta\": {\n            \"num_b\": len(le_b.classes_),\n            \"num_s\": len(le_s.classes_),\n            \"num_u\": len(le_u.classes_),\n            \"feats\": feat_cols,\n        },\n    }\n\n# =========================================================\n# 2. MODEL (UNCHANGED CORE STRUCTURE)\n# =========================================================\nclass PIGD_Block(nn.Module):\n    def __init__(self, d, c):\n        super().__init__()\n        self.ln = nn.LayerNorm(d)\n        self.g = nn.Linear(c, d)\n        self.b = nn.Linear(c, d)\n        self.gate = nn.Sequential(nn.Linear(d, d), nn.Sigmoid())\n        self.net = nn.Sequential(\n            nn.Linear(d, d * 2),\n            nn.SiLU(),\n            nn.Linear(d * 2, d)\n        )\n\n    def forward(self, x, ctx):\n        h = self.ln(x)\n        h = h * (1 + self.g(ctx)) + self.b(ctx)\n        return x + self.gate(h) * self.net(h)\n\nclass PIGD_Model(nn.Module):\n    def __init__(self, inp, nb, ns, nu):\n        super().__init__()\n        self.b = nn.Embedding(nb, 64)\n        self.s = nn.Embedding(ns, 16)\n        self.u = nn.Embedding(nu, 16)\n        self.t = nn.Sequential(nn.Linear(1, 32), nn.SiLU())\n\n        ctx_dim = 64 + 16 + 16 + 32\n        self.inp = nn.Linear(inp, 256)\n        self.blocks = nn.ModuleList([PIGD_Block(256, ctx_dim) for _ in range(6)])\n        self.out = nn.Linear(256, inp)\n\n    def forward(self, x, t, b, s, u):\n        te = self.t(t.float().unsqueeze(1) / 100)\n        ctx = torch.cat([self.b(b), self.s(s), self.u(u), te], dim=1)\n        h = self.inp(x)\n        for blk in self.blocks:\n            h = blk(h, ctx)\n        return self.out(h)\n\n# =========================================================\n# SAFE SAVITZKY–GOLAY\n# =========================================================\ndef safe_savgol(x):\n    x = np.asarray(x, dtype=np.float32)\n    if len(x) < 7 or np.allclose(x, x[0]):\n        return x\n    w = min(15, len(x)//2*2 + 1)\n    try:\n        return savgol_filter(x, w, 2, mode=\"interp\")\n    except:\n        return x\n\n# =========================================================\n# 3. EXECUTION WITH K-FOLD = 3 + EPOCHS\n# =========================================================\ndef run_piml_engine(data):\n    meta = data[\"meta\"]\n    x_all, b_all, s_all, u_all = data[\"train\"]\n\n    kf = KFold(n_splits=3, shuffle=True, random_state=42)\n\n    aucs, prs = [], []\n\n    EPOCHS = 25\n    LR = 7e-4\n\n    for fold, (idx, _) in enumerate(kf.split(x_all)):\n        print(f\"\\n🔥 Fold {fold+1}/3 Training...\")\n\n        model = PIGD_Model(\n            len(meta[\"feats\"]),\n            meta[\"num_b\"],\n            meta[\"num_s\"],\n            meta[\"num_u\"]\n        ).to(device)\n\n        opt = optim.AdamW(model.parameters(), lr=LR)\n        scaler = torch.cuda.amp.GradScaler()\n\n        # Softer diffusion noise\n        alphas = torch.cumprod(\n            1 - torch.linspace(5e-5, 0.01, 100, device=device), 0\n        )\n\n        loader = DataLoader(\n            TensorDataset(x_all[idx], b_all[idx], s_all[idx], u_all[idx]),\n            batch_size=2048,\n            shuffle=True\n        )\n\n        model.train()\n        for epoch in range(EPOCHS):\n            for x0, b0, s0, u0 in loader:\n                x0, b0, s0, u0 = (\n                    x0.to(device),\n                    b0.to(device),\n                    s0.to(device),\n                    u0.to(device),\n                )\n\n                t = torch.randint(0, 100, (x0.size(0),), device=device)\n                a = alphas[t].unsqueeze(1)\n                xt = torch.sqrt(a) * x0 + torch.sqrt(1 - a) * torch.randn_like(x0)\n\n                opt.zero_grad()\n                with torch.cuda.amp.autocast():\n                    loss = nn.MSELoss()(model(xt, t, b0, s0, u0), x0)\n\n                scaler.scale(loss).backward()\n                scaler.step(opt)\n                scaler.update()\n\n            print(f\"Epoch {epoch+1}/{EPOCHS} | Loss {loss.item():.4f}\")\n\n        # Inference\n        model.eval()\n        x, b, s, u, y, bid = data[\"test\"]\n        x, b, s, u = x.to(device), b.to(device), s.to(device), u.to(device)\n\n        with torch.no_grad():\n            t = torch.full((x.size(0),), 25, device=device)\n            a = alphas[t].unsqueeze(1)\n            xt = torch.sqrt(a) * x + torch.sqrt(1 - a) * torch.randn_like(x)\n            r = model(xt, t, b, s, u)\n\n            shock_err   = torch.abs(x[:, 0] - r[:, 0])\n            kinetic_err = torch.abs(x[:, 2] - r[:, 2])\n            thermal_err = torch.abs(x[:, 1] - r[:, 1])\n\n            err = (shock_err * (1 + 0.5 * kinetic_err) + 0.3 * thermal_err).cpu().numpy()\n\n        df = pd.DataFrame({\"e\": err, \"b\": bid})\n        smooth = df.groupby(\"b\")[\"e\"].transform(safe_savgol)\n        smooth = np.nan_to_num(smooth, 0)\n\n        auc = roc_auc_score(y, smooth)\n        pr = average_precision_score(y, smooth)\n\n        print(f\"Fold {fold+1} ROC-AUC: {auc:.4f}\")\n        print(f\"Fold {fold+1} PR-AUC : {pr:.4f}\")\n\n        aucs.append(auc)\n        prs.append(pr)\n\n    print(\"\\n🏆 FINAL K-FOLD (k=3) RESULTS\")\n    print(f\"Mean ROC-AUC: {np.mean(aucs):.4f}\")\n    print(f\"Mean PR-AUC : {np.mean(prs):.4f}\")\n\n# =========================================================\n# RUN (KAGGLE PATH)\n# =========================================================\nif __name__ == \"__main__\":\n    path = \"/kaggle/input/energy-anomaly-detection/train_features.csv\"\n    if os.path.exists(path):\n        data = prepare_physics_sota_data(path)\n        run_piml_engine(data)\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-28T18:37:40.557396Z","iopub.execute_input":"2026-01-28T18:37:40.55782Z","iopub.status.idle":"2026-01-28T18:57:07.283546Z","shell.execute_reply.started":"2026-01-28T18:37:40.557786Z","shell.execute_reply":"2026-01-28T18:57:07.282696Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# =========================================================\n# Kaggle Environment & Imports\n# =========================================================\nimport os\nimport warnings\nwarnings.filterwarnings(\"ignore\")\n\nimport numpy as np\nimport pandas as pd\n\nimport torch\nfrom torch import nn, optim\nfrom torch.utils.data import DataLoader, TensorDataset\n\nfrom sklearn.preprocessing import LabelEncoder\nfrom sklearn.metrics import roc_auc_score, average_precision_score\nfrom sklearn.model_selection import KFold\n\nfrom scipy.signal import savgol_filter\n\n# =========================================================\n# DEVICE\n# =========================================================\ndevice = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\nprint(\"Device:\", device)\nif device.type == \"cuda\":\n    print(\"GPU:\", torch.cuda.get_device_name(0))\n\n# =========================================================\n# 1. DATA PREPARATION (ORIGINAL LOGIC PRESERVED)\n# =========================================================\ndef prepare_physics_sota_data(path):\n    print(\"🧬 Engineering Thermal Inertia & Kinetic Momentum Manifolds...\")\n\n    df = pd.read_csv(path)\n\n    # Fixed timestamp parsing\n    df[\"timestamp\"] = pd.to_datetime(df[\"timestamp\"], format=\"mixed\", errors=\"coerce\")\n    df = df.dropna(subset=[\"timestamp\"])\n    df = df.sort_values([\"building_id\", \"timestamp\"])\n\n    # Thermodynamic baseline\n    df[\"meter_log\"] = np.log1p(df[\"meter_reading\"])\n    df[\"hour\"] = df[\"timestamp\"].dt.hour\n    df[\"dow\"] = df[\"timestamp\"].dt.dayofweek\n\n    norm = (\n        df[df[\"anomaly\"] == 0]\n        .groupby([\"building_id\", \"dow\", \"hour\"])[\"meter_log\"]\n        .agg([\"median\", \"std\"])\n        .reset_index()\n    )\n    norm.columns = [\"building_id\", \"dow\", \"hour\", \"h_med\", \"h_std\"]\n    df = df.merge(norm, on=[\"building_id\", \"dow\", \"hour\"], how=\"left\")\n\n    df[\"h_med\"] = df[\"h_med\"].fillna(\n        df.groupby(\"building_id\")[\"meter_log\"].transform(\"median\")\n    )\n    df[\"h_std\"] = df[\"h_std\"].replace(0, 1).fillna(1)\n\n    # Physics features\n    df[\"temp_inertia\"] = (\n        df.groupby(\"building_id\")[\"air_temperature\"]\n        .transform(lambda x: x.rolling(6, min_periods=1).mean())\n    )\n    df[\"thermal_p\"] = df[\"air_temperature\"] - df[\"temp_inertia\"]\n\n    df[\"shock\"] = df[\"meter_log\"] - df[\"h_med\"]\n    df[\"shock_v\"] = df.groupby(\"building_id\")[\"shock\"].diff().fillna(0)\n    df[\"kinetic_s\"] = df.groupby(\"building_id\")[\"shock_v\"].diff().fillna(0)\n\n    # Robust scaling\n    df[\"shock_scaled\"] = df.groupby(\"building_id\")[\"shock\"].transform(\n        lambda x: (x - x.median()) /\n        (x.quantile(0.75) - x.quantile(0.25) + 1e-6)\n    )\n\n    feat_cols = [\"shock_scaled\", \"thermal_p\", \"kinetic_s\", \"air_temperature\"]\n\n    df[feat_cols] = df[feat_cols].replace([np.inf, -np.inf], 0).fillna(0)\n    df[feat_cols] = df[feat_cols].clip(-25, 25)\n\n    # Encoding\n    le_b, le_s, le_u = LabelEncoder(), LabelEncoder(), LabelEncoder()\n    df[\"b_idx\"] = le_b.fit_transform(df[\"building_id\"])\n    df[\"s_idx\"] = le_s.fit_transform(df[\"site_id\"])\n    df[\"u_idx\"] = le_u.fit_transform(df[\"primary_use\"])\n\n    train_df = df[df[\"anomaly\"] == 0].sample(frac=0.45, random_state=42)\n    test_df = df.drop(train_df.index).sort_values([\"building_id\", \"timestamp\"])\n\n    return {\n        \"train\": (\n            torch.tensor(train_df[feat_cols].values, dtype=torch.float32),\n            torch.tensor(train_df[\"b_idx\"].values),\n            torch.tensor(train_df[\"s_idx\"].values),\n            torch.tensor(train_df[\"u_idx\"].values),\n        ),\n        \"test\": (\n            torch.tensor(test_df[feat_cols].values, dtype=torch.float32),\n            torch.tensor(test_df[\"b_idx\"].values),\n            torch.tensor(test_df[\"s_idx\"].values),\n            torch.tensor(test_df[\"u_idx\"].values),\n            test_df[\"anomaly\"].values.astype(np.int32),\n            test_df[\"building_id\"].values,\n        ),\n        \"meta\": {\n            \"num_b\": len(le_b.classes_),\n            \"num_s\": len(le_s.classes_),\n            \"num_u\": len(le_u.classes_),\n            \"feats\": feat_cols,\n        },\n    }\n\n# =========================================================\n# 2. MODEL (UNCHANGED CORE STRUCTURE)\n# =========================================================\nclass PIGD_Block(nn.Module):\n    def __init__(self, d, c):\n        super().__init__()\n        self.ln = nn.LayerNorm(d)\n        self.g = nn.Linear(c, d)\n        self.b = nn.Linear(c, d)\n        self.gate = nn.Sequential(nn.Linear(d, d), nn.Sigmoid())\n        self.net = nn.Sequential(\n            nn.Linear(d, d * 2),\n            nn.SiLU(),\n            nn.Linear(d * 2, d)\n        )\n\n    def forward(self, x, ctx):\n        h = self.ln(x)\n        h = h * (1 + self.g(ctx)) + self.b(ctx)\n        return x + self.gate(h) * self.net(h)\n\nclass PIGD_Model(nn.Module):\n    def __init__(self, inp, nb, ns, nu):\n        super().__init__()\n        self.b = nn.Embedding(nb, 64)\n        self.s = nn.Embedding(ns, 16)\n        self.u = nn.Embedding(nu, 16)\n        self.t = nn.Sequential(nn.Linear(1, 32), nn.SiLU())\n\n        ctx_dim = 64 + 16 + 16 + 32\n        self.inp = nn.Linear(inp, 256)\n        self.blocks = nn.ModuleList([PIGD_Block(256, ctx_dim) for _ in range(6)])\n        self.out = nn.Linear(256, inp)\n\n    def forward(self, x, t, b, s, u):\n        te = self.t(t.float().unsqueeze(1) / 100)\n        ctx = torch.cat([self.b(b), self.s(s), self.u(u), te], dim=1)\n        h = self.inp(x)\n        for blk in self.blocks:\n            h = blk(h, ctx)\n        return self.out(h)\n\n# =========================================================\n# SAFE SAVITZKY–GOLAY\n# =========================================================\ndef safe_savgol(x):\n    x = np.asarray(x, dtype=np.float32)\n    if len(x) < 7 or np.allclose(x, x[0]):\n        return x\n    w = min(15, len(x)//2*2 + 1)\n    try:\n        return savgol_filter(x, w, 2, mode=\"interp\")\n    except:\n        return x\n\n# =========================================================\n# 3. EXECUTION WITH K-FOLD = 3 + EPOCHS\n# =========================================================\ndef run_piml_engine(data):\n    meta = data[\"meta\"]\n    x_all, b_all, s_all, u_all = data[\"train\"]\n\n    kf = KFold(n_splits=3, shuffle=True, random_state=42)\n\n    aucs, prs = [], []\n\n    EPOCHS = 25\n    LR = 7e-4\n\n    for fold, (idx, _) in enumerate(kf.split(x_all)):\n        print(f\"\\n🔥 Fold {fold+1}/3 Training...\")\n\n        model = PIGD_Model(\n            len(meta[\"feats\"]),\n            meta[\"num_b\"],\n            meta[\"num_s\"],\n            meta[\"num_u\"]\n        ).to(device)\n\n        opt = optim.AdamW(model.parameters(), lr=LR)\n        scaler = torch.cuda.amp.GradScaler()\n\n        # Softer diffusion noise\n        alphas = torch.cumprod(\n            1 - torch.linspace(5e-5, 0.01, 100, device=device), 0\n        )\n\n        loader = DataLoader(\n            TensorDataset(x_all[idx], b_all[idx], s_all[idx], u_all[idx]),\n            batch_size=2048,\n            shuffle=True\n        )\n\n        model.train()\n        for epoch in range(EPOCHS):\n            for x0, b0, s0, u0 in loader:\n                x0, b0, s0, u0 = (\n                    x0.to(device),\n                    b0.to(device),\n                    s0.to(device),\n                    u0.to(device),\n                )\n\n                t = torch.randint(0, 100, (x0.size(0),), device=device)\n                a = alphas[t].unsqueeze(1)\n                xt = torch.sqrt(a) * x0 + torch.sqrt(1 - a) * torch.randn_like(x0)\n\n                opt.zero_grad()\n                with torch.cuda.amp.autocast():\n                    loss = nn.MSELoss()(model(xt, t, b0, s0, u0), x0)\n\n                scaler.scale(loss).backward()\n                scaler.step(opt)\n                scaler.update()\n\n            print(f\"Epoch {epoch+1}/{EPOCHS} | Loss {loss.item():.4f}\")\n\n        # Inference\n        model.eval()\n        x, b, s, u, y, bid = data[\"test\"]\n        x, b, s, u = x.to(device), b.to(device), s.to(device), u.to(device)\n\n        with torch.no_grad():\n            t = torch.full((x.size(0),), 25, device=device)\n            a = alphas[t].unsqueeze(1)\n            xt = torch.sqrt(a) * x + torch.sqrt(1 - a) * torch.randn_like(x)\n            r = model(xt, t, b, s, u)\n\n            shock_err   = torch.abs(x[:, 0] - r[:, 0])\n            kinetic_err = torch.abs(x[:, 2] - r[:, 2])\n            thermal_err = torch.abs(x[:, 1] - r[:, 1])\n\n            err = (shock_err * (1 + 0.5 * kinetic_err) + 0.3 * thermal_err).cpu().numpy()\n\n        df = pd.DataFrame({\"e\": err, \"b\": bid})\n        smooth = df.groupby(\"b\")[\"e\"].transform(safe_savgol)\n        smooth = np.nan_to_num(smooth, 0)\n\n        auc = roc_auc_score(y, smooth)\n        pr = average_precision_score(y, smooth)\n\n        print(f\"Fold {fold+1} ROC-AUC: {auc:.4f}\")\n        print(f\"Fold {fold+1} PR-AUC : {pr:.4f}\")\n\n        aucs.append(auc)\n        prs.append(pr)\n\n    print(\"\\n🏆 FINAL K-FOLD (k=3) RESULTS\")\n    print(f\"Mean ROC-AUC: {np.mean(aucs):.4f}\")\n    print(f\"Mean PR-AUC : {np.mean(prs):.4f}\")\n\n# =========================================================\n# RUN (KAGGLE PATH)\n# =========================================================\nif __name__ == \"__main__\":\n    path = \"/kaggle/input/energy-anomaly-detection/train_features.csv\"\n    if os.path.exists(path):\n        data = prepare_physics_sota_data(path)\n        run_piml_engine(data)\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-28T19:10:04.380244Z","iopub.execute_input":"2026-01-28T19:10:04.380585Z","iopub.status.idle":"2026-01-28T19:28:40.187043Z","shell.execute_reply.started":"2026-01-28T19:10:04.380558Z","shell.execute_reply":"2026-01-28T19:28:40.186261Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# =========================================================\n# 0. IMPORTS\n# =========================================================\nimport os\nimport warnings\nwarnings.filterwarnings(\"ignore\")\n\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\n\nimport torch\nfrom torch import nn, optim\nfrom torch.utils.data import DataLoader, TensorDataset\n\nfrom sklearn.metrics import (\n    roc_auc_score,\n    average_precision_score,\n    precision_score,\n    recall_score,\n    roc_curve,\n    precision_recall_curve\n)\n\nfrom scipy.signal import savgol_filter\n\ndevice = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\nprint(\"Device:\", device)\ntorch.cuda.empty_cache()\n\n# =========================================================\n# PARAMETERS\n# =========================================================\nDIFF_STEPS = 100\nINFER_T = 25\nMC_PASSES = 7\n\n# =========================================================\n# 1. DATA PREPARATION\n# =========================================================\ndef prepare_data(path):\n\n    df = pd.read_csv(path)\n\n    df[\"timestamp\"] = pd.to_datetime(df[\"timestamp\"], format=\"mixed\", errors=\"coerce\")\n    df = df.dropna(subset=[\"timestamp\"])\n    df = df.sort_values([\"building_id\", \"timestamp\"])\n\n    df[\"meter_log\"] = np.log1p(df[\"meter_reading\"])\n    df[\"hour\"] = df[\"timestamp\"].dt.hour\n    df[\"dow\"] = df[\"timestamp\"].dt.dayofweek\n\n    norm = (\n        df[df[\"anomaly\"] == 0]\n        .groupby([\"building_id\",\"dow\",\"hour\"])[\"meter_log\"]\n        .agg([\"median\",\"std\"])\n        .reset_index()\n    )\n\n    norm.columns = [\"building_id\",\"dow\",\"hour\",\"h_med\",\"h_std\"]\n    df = df.merge(norm, on=[\"building_id\",\"dow\",\"hour\"], how=\"left\")\n\n    df[\"h_med\"] = df[\"h_med\"].fillna(\n        df.groupby(\"building_id\")[\"meter_log\"].transform(\"median\")\n    )\n\n    df[\"temp_inertia\"] = (\n        df.groupby(\"building_id\")[\"air_temperature\"]\n        .transform(lambda x: x.rolling(6,min_periods=1).mean())\n    )\n\n    df[\"thermal_p\"] = df[\"air_temperature\"] - df[\"temp_inertia\"]\n\n    df[\"shock\"] = df[\"meter_log\"] - df[\"h_med\"]\n    df[\"shock_v\"] = df.groupby(\"building_id\")[\"shock\"].diff().fillna(0)\n    df[\"kinetic_s\"] = df.groupby(\"building_id\")[\"shock_v\"].diff().fillna(0)\n\n    df[\"shock_scaled\"] = df.groupby(\"building_id\")[\"shock\"].transform(\n        lambda x: (x - x.median()) /\n        (x.quantile(0.75) - x.quantile(0.25) + 1e-6)\n    )\n\n    feat_cols = [\"shock_scaled\",\"thermal_p\",\"kinetic_s\",\"air_temperature\"]\n    df[feat_cols] = df[feat_cols].replace([np.inf,-np.inf],0).fillna(0)\n\n    return df, feat_cols\n\n# =========================================================\n# 2. MODEL (STRONG BASELINE)\n# =========================================================\nclass PIGD(nn.Module):\n    def __init__(self, input_dim):\n        super().__init__()\n        self.net = nn.Sequential(\n            nn.Linear(input_dim,128),\n            nn.SiLU(),\n            nn.Linear(128,128),\n            nn.SiLU(),\n            nn.Linear(128,input_dim)\n        )\n\n    def forward(self,x):\n        return torch.nan_to_num(self.net(x))\n\n# =========================================================\n# 3. SAFE SAVITZKY-GOLAY\n# =========================================================\ndef safe_savgol(x):\n    x = np.asarray(x, dtype=np.float64)\n\n    if len(x) < 7:\n        return x\n\n    w = min(11, len(x))\n    if w % 2 == 0:\n        w -= 1\n\n    if w < 3:\n        return x\n\n    try:\n        return savgol_filter(x, window_length=w, polyorder=2, mode=\"interp\")\n    except:\n        return x\n\n# =========================================================\n# 4. TRAIN + MONTE CARLO DIFFUSION\n# =========================================================\ndef train_model(df, feat_cols):\n\n    X = torch.tensor(df[feat_cols].values, dtype=torch.float32)\n    y = df[\"anomaly\"].values\n    b_ids = df[\"building_id\"].values\n\n    normal_mask = y == 0\n    X_train = X[normal_mask]\n\n    loader = DataLoader(\n        TensorDataset(X_train),\n        batch_size=1024,\n        shuffle=True\n    )\n\n    model = PIGD(len(feat_cols)).to(device)\n    opt = optim.AdamW(model.parameters(), lr=1e-3)\n\n    alphas = torch.linspace(1e-4,0.02,DIFF_STEPS).to(device)\n\n    print(\"🔥 Training Strong PIGD\")\n\n    for epoch in range(20):\n        model.train()\n        for (xb,) in loader:\n            xb = xb.to(device)\n\n            t = torch.randint(0,DIFF_STEPS,(xb.size(0),),device=device)\n            noise = torch.randn_like(xb)\n            alpha = alphas[t].unsqueeze(1)\n\n            xt = torch.sqrt(1-alpha)*xb + torch.sqrt(alpha)*noise\n\n            pred = model(xt)\n            loss = nn.MSELoss()(pred, xb)\n\n            opt.zero_grad()\n            loss.backward()\n            torch.nn.utils.clip_grad_norm_(model.parameters(),1.0)\n            opt.step()\n\n        print(f\"Epoch {epoch+1} | Loss {loss.item():.6f}\")\n\n    # =============================\n    # MONTE CARLO INFERENCE\n    # =============================\n    model.eval()\n\n    full_loader = DataLoader(\n        TensorDataset(X),\n        batch_size=2048,\n        shuffle=False\n    )\n\n    scores = np.zeros(len(X))\n\n    with torch.no_grad():\n        idx = 0\n        for (xb,) in full_loader:\n            xb = xb.to(device)\n            batch_scores = torch.zeros(xb.size(0),device=device)\n\n            for _ in range(MC_PASSES):\n                t = torch.full((xb.size(0),),INFER_T,device=device)\n                noise = torch.randn_like(xb)\n                alpha = alphas[t].unsqueeze(1)\n\n                xt = torch.sqrt(1-alpha)*xb + torch.sqrt(alpha)*noise\n                pred = model(xt)\n\n                err = (\n                    (pred[:,0]-xb[:,0])**2 +\n                    0.8*(pred[:,2]-xb[:,2])**2 +\n                    0.5*(pred[:,1]-xb[:,1])**2\n                )\n\n                batch_scores += err\n\n            batch_scores /= MC_PASSES\n\n            scores[idx:idx+xb.size(0)] = batch_scores.cpu().numpy()\n            idx += xb.size(0)\n\n    # =============================\n    # SAFE PER-BUILDING SMOOTHING\n    # =============================\n    df_temp = pd.DataFrame({\n        \"score\": scores,\n        \"building\": b_ids\n    })\n\n    df_temp[\"score\"] = df_temp.groupby(\"building\")[\"score\"].transform(safe_savgol)\n\n    scores = np.nan_to_num(df_temp[\"score\"].values)\n    scores = np.log1p(scores)\n\n    return scores, y, b_ids\n\n# =========================================================\n# 5. EVALUATION\n# =========================================================\ndef evaluate(scores, y, b_ids):\n\n    auc = roc_auc_score(y,scores)\n    pr = average_precision_score(y,scores)\n\n    thresh = np.percentile(scores,95)\n    preds = (scores>thresh).astype(int)\n\n    precision = precision_score(y,preds)\n    recall = recall_score(y,preds)\n    fpr = np.sum((preds==1)&(y==0)) / np.sum(y==0)\n\n    print(\"\\n🏆 FINAL RESULTS\")\n    print(\"ROC-AUC :\",auc)\n    print(\"PR-AUC  :\",pr)\n    print(\"Precision:\",precision)\n    print(\"Recall:\",recall)\n    print(\"TPR:\",recall)\n    print(\"FPR:\",fpr)\n\n    fpr_c,tpr_c,_ = roc_curve(y,scores)\n    plt.plot(fpr_c,tpr_c)\n    plt.title(\"ROC Curve\")\n    plt.show()\n\n    p_c,r_c,_ = precision_recall_curve(y,scores)\n    plt.plot(r_c,p_c)\n    plt.title(\"PR Curve\")\n    plt.show()\n\n# =========================================================\n# RUN\n# =========================================================\npath = \"/kaggle/input/energy-anomaly-detection/train_features.csv\"\n\ndf, feats = prepare_data(path)\nscores, y, b_ids = train_model(df, feats)\nevaluate(scores, y, b_ids)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-03-03T23:06:28.273771Z","iopub.execute_input":"2026-03-03T23:06:28.274509Z","iopub.status.idle":"2026-03-03T23:10:04.307043Z","shell.execute_reply.started":"2026-03-03T23:06:28.274471Z","shell.execute_reply":"2026-03-03T23:10:04.305752Z"}},"outputs":[],"execution_count":null}]}