Added full Markov-model pipeline (#11)

* added full markov model

* comments

* use of config

* test loader

* added tests

Author: PhilippSchmelter

CHANGELOG.md

@@ -8,6 +8,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ### Added
 - `CODEOWNERS` file and dependabot automation [#7](https://github.com/ie3-institute/simonaMarkovLoad/issues/7)
+- Added full Markov-model pipeline [#10](https://github.com/ie3-institute/simonaMarkovLoad/issues/10)
 
 ### Changed
 - Compute instantaneous kW from cumulative kWh via 15-minute differencing [#1](https://github.com/ie3-institute/simonaMarkovLoad/issues/1)

config.yml

@@ -3,3 +3,7 @@ input:
   timestamp_col: "Zeitstempel"
   value_col: "Messwert"
   factor: 4         # 15-min → kW depends on sampling
+
+model:
+  n_states: 10
+  laplace_alpha: 1.0

src/config.py

@@ -1,4 +1,3 @@
-# src/config.py
 import os
 from pathlib import Path
 

src/main.py

@@ -1,30 +1,31 @@
 import matplotlib.pyplot as plt
+import numpy as np
 
+from src.markov.transition_counts import build_transition_counts
+from src.markov.transitions import build_transition_matrices
 from src.preprocessing.loader import load_timeseries
 
+df = load_timeseries(normalize=True, discretize=True)
 
-def plot_state_distribution(df):
 
-    counts = df["state"].value_counts().sort_index()
+counts = build_transition_counts(df)
+probs = build_transition_matrices(df, alpha=1.0)
 
-    plt.figure()
-    plt.bar(counts.index, counts.values)
-    plt.xlabel("State")
-    plt.ylabel("Anzahl Einträge")
-    plt.title("Verteilung der Einträge nach State")
-    plt.xticks(counts.index)
-    plt.show()
+print("counts shape :", counts.shape)
+print("probs  shape :", probs.shape)
 
 
-def main():
-    df = load_timeseries()
-    print(df)
-    df_norm = load_timeseries(normalize=True)
-    print(df_norm)
-    df_disc = load_timeseries(normalize=True, discretize=True)
-    print(df_disc)
-    plot_state_distribution(df_disc)
+active_buckets = np.where(counts.sum(axis=(1, 2)) > 0)[0]
+bucket = int(active_buckets[0]) if active_buckets.size else 0
+print(f"\nUsing bucket {bucket}")
 
 
-if __name__ == "__main__":
-    main()
+print("row sums :", probs[bucket].sum(axis=1))
+
+plt.imshow(probs[bucket], aspect="auto")
+plt.title(f"Bucket {bucket} – transition probabilities")
+plt.xlabel("state t+1")
+plt.ylabel("state t")
+plt.colorbar()
+plt.tight_layout()
+plt.show()

src/markov/__init__.py

@@ -0,0 +1,10 @@
+from .buckets import NUM_BUCKETS, bucket_id
+from .transition_counts import build_transition_counts
+from .transitions import build_transition_matrices
+
+__all__ = [
+    "bucket_id",
+    "NUM_BUCKETS",
+    "build_transition_counts",
+    "build_transition_matrices",
+]

src/markov/_core.py

@@ -0,0 +1,21 @@
+import numpy as np
+import pandas as pd
+
+from src.config import CONFIG
+
+from .buckets import NUM_BUCKETS
+
+N_STATES = int(CONFIG["model"]["n_states"])
+
+
+def _transition_counts(
+    df: pd.DataFrame, *, state_col="state", bucket_col="bucket", dtype=np.uint32
+) -> np.ndarray:
+    df = df.sort_values("timestamp")
+    s_t = df[state_col].to_numpy(dtype=int)[:-1]
+    s_tp1 = df[state_col].to_numpy(dtype=int)[1:]
+    buckets = df[bucket_col].to_numpy(dtype=int)[:-1]
+
+    counts = np.zeros((NUM_BUCKETS, N_STATES, N_STATES), dtype=dtype)
+    np.add.at(counts, (buckets, s_t, s_tp1), 1)
+    return counts

src/markov/buckets.py

@@ -0,0 +1,44 @@
+import pandas as pd
+
+_MONTH_FACTOR = 96 * 2
+_WEEKEND_FACTOR = 96
+_NUM_MONTHS = 12
+_NUM_QH = 96
+NUM_BUCKETS = _NUM_MONTHS * 2 * _NUM_QH  # 2 304
+
+
+def _is_weekend(ts):
+    if isinstance(ts, pd.Series):
+        return ts.dt.dayofweek >= 5
+    if isinstance(ts, pd.DatetimeIndex):
+        return ts.dayofweek >= 5
+    return ts.dayofweek >= 5
+
+
+def bucket_id(ts):
+    if isinstance(ts, pd.Series):
+        weekend = _is_weekend(ts).astype(int)
+        qh = ts.dt.hour * 4 + ts.dt.minute // 15
+        month = ts.dt.month - 1
+    elif isinstance(ts, pd.DatetimeIndex):
+        weekend = _is_weekend(ts).astype(int)
+        qh = ts.hour * 4 + ts.minute // 15
+        month = ts.month - 1
+    else:
+        weekend = int(_is_weekend(ts))
+        qh = ts.hour * 4 + ts.minute // 15
+        month = ts.month - 1
+
+    return month * _MONTH_FACTOR + weekend * _WEEKEND_FACTOR + qh
+
+
+def assign_buckets(
+    df: pd.DataFrame,
+    *,
+    ts_col: str = "timestamp",
+    bucket_col: str = "bucket",
+    inplace: bool = False,
+) -> pd.DataFrame:
+    tgt = df if inplace else df.copy()
+    tgt[bucket_col] = bucket_id(tgt[ts_col]).astype("uint16")
+    return tgt

src/markov/model.py

@@ -0,0 +1,31 @@
+import numpy as np
+import pandas as pd
+
+from src.config import CONFIG
+
+from .buckets import NUM_BUCKETS
+
+N_STATES = int(CONFIG["model"]["n_states"])
+
+
+def build_transition_counts(
+    df: pd.DataFrame,
+    *,
+    state_col: str = "state",
+    bucket_col: str = "bucket",
+    dtype=np.uint32,
+) -> np.ndarray:
+    """
+    Absolute transition counts:
+        C[b, i, j] = # of times state_t=i  → state_{t+1}=j in bucket b
+    Shape = (2 304, 10, 10).
+    """
+    df = df.sort_values("timestamp")
+
+    s_t = df[state_col].to_numpy(dtype=int)[:-1]
+    s_tp1 = df[state_col].to_numpy(dtype=int)[1:]
+    buckets = df[bucket_col].to_numpy(dtype=int)[:-1]
+
+    counts = np.zeros((NUM_BUCKETS, N_STATES, N_STATES), dtype=dtype)
+    np.add.at(counts, (buckets, s_t, s_tp1), 1)
+    return counts

src/markov/transition_counts.py

@@ -0,0 +1,15 @@
+import numpy as np
+import pandas as pd
+
+from src.config import CONFIG
+
+from ._core import _transition_counts
+
+alpha = float(CONFIG["model"]["laplace_alpha"])
+
+
+def build_transition_matrices(df: pd.DataFrame, *, dtype=np.float32) -> np.ndarray:
+    counts = _transition_counts(df, dtype=dtype)
+    counts += alpha
+    counts /= counts.sum(axis=2, keepdims=True)
+    return counts.astype(dtype)