From 7f11dce9b19223e0b91781fdbe6cf69cdd4a89f9 Mon Sep 17 00:00:00 2001 From: Claude Date: Thu, 16 Jul 2026 20:58:45 +0000 Subject: [PATCH 1/8] Start batch modernization: offline batch PCA, dataset loaders, converter completion Opens the Phase 1 work described in the batch data processing modernization plan: an offline batchwise (BWU) PCA model with batch-level diagnostics, dataset loader functions for the bundled nylon/dryer/fake batch data, completion of the three unimplemented data_input converters, and a sweep of small lint/dead-code cleanups. Co-Authored-By: Claude Fable 5 Claude-Session: https://claude.ai/code/session_01SyEpexCyqHn1rwTMYQSiP6 From bbba73980d8b9ad2e230f242c17d92e051eda6a6 Mon Sep 17 00:00:00 2001 From: Claude Date: Thu, 16 Jul 2026 21:03:08 +0000 Subject: [PATCH 2/8] Add loader functions for the bundled batch datasets load_nylon(), load_dryer(), and load_batch_fake_data() return the standard batch-data dictionary directly, following the loader-function precedent in experiments/datasets.py. Until now every consumer had to locate the CSVs by path and convert them by hand; the book's worked examples need a one-line reproducible load. Co-Authored-By: Claude Fable 5 Claude-Session: https://claude.ai/code/session_01SyEpexCyqHn1rwTMYQSiP6 --- src/process_improve/batch/datasets.py | 143 ++++++++++++++++++++++++++ tests/batch/test_datasets.py | 45 ++++++++ 2 files changed, 188 insertions(+) create mode 100644 src/process_improve/batch/datasets.py create mode 100644 tests/batch/test_datasets.py diff --git a/src/process_improve/batch/datasets.py b/src/process_improve/batch/datasets.py new file mode 100644 index 00000000..5d9d702a --- /dev/null +++ b/src/process_improve/batch/datasets.py @@ -0,0 +1,143 @@ +# (c) Kevin Dunn, 2010-2026. MIT License. Based on own private work over the years. +"""Loader functions for the batch datasets bundled with the package. + +Each loader returns the standard batch-data dictionary used throughout +:mod:`process_improve.batch`: keys are batch identifiers, values are per-batch +dataframes with identical, all-numeric columns (one column per tag). See +:mod:`process_improve.batch.data_input` for the format definitions and +converters to the melted and wide representations. +""" + +from __future__ import annotations + +from pathlib import Path +from typing import TYPE_CHECKING + +import pandas as pd + +from .data_input import melted_to_dict + +if TYPE_CHECKING: + from collections.abc import Hashable + +_DATASETS_DIR = Path(__file__).resolve().parents[1] / "datasets" / "batch" + + +def _load_melted_csv( + filename: str, batch_id_col: str, drop_columns: list[str] | None = None +) -> dict[Hashable, pd.DataFrame]: + """Read a bundled melted CSV and return it as a batch-data dictionary. + + The ``batch_id_col`` column is consumed as the dictionary key and dropped + from each per-batch dataframe (the key carries that information), along + with any extra ``drop_columns``. Row indexes are reset per batch, so each + batch counts its own samples from zero. + """ + melted = pd.read_csv(_DATASETS_DIR / filename) + batches = melted_to_dict(melted, batch_id_col=batch_id_col) + to_drop = [batch_id_col, *(drop_columns or [])] + return {batch_id: batch.drop(columns=to_drop).reset_index(drop=True) for batch_id, batch in batches.items()} + + +def load_nylon() -> dict[Hashable, pd.DataFrame]: + """Return the nylon autoclave reactor batch dataset. + + Trajectory data from an industrial nylon polymerization autoclave, + used widely in the batch analysis and monitoring literature. Variables + ``Tag01`` to ``Tag10`` are temperatures, pressures, and flows recorded + during each batch. Batch durations vary slightly (113 to 135 samples), + so resample or align the batches to a common length before unfolding + (see :func:`process_improve.batch.resample_to_reference`). + + Returns + ------- + dict[Hashable, pd.DataFrame] + Standard batch-data dictionary: 57 batches, each a dataframe of + 10 numeric tag columns. + + Source + ------ + Kassidas, A., "Fault Detection and Diagnosis in Dynamic Multivariable + Chemical Processes Using Speech Recognition Methods", PhD thesis, + McMaster University, 1997. Also analyzed in Wold, Kettaneh-Wold, + MacGregor and Dunn, "Batch Process Modeling and MSPC", Comprehensive + Chemometrics, Elsevier, 2009. + + Examples + -------- + >>> from process_improve.batch.datasets import load_nylon + >>> batches = load_nylon() + >>> len(batches) + 57 + """ + return _load_melted_csv("nylon.csv", batch_id_col="batch_id") + + +def load_dryer() -> dict[Hashable, pd.DataFrame]: + """Return the batch dryer dataset. + + Trajectory data from an industrial batch drying process (a critical step + in the manufacture of an agricultural chemical). Each batch records ten + process tags plus ``ClockTime``, the wall-time sample counter: + + - ``CollectorTankLevel``: level of the solvent collector tank + - ``DifferentialPressure``: differential pressure in the dryer + - ``DryerPressure``: pressure in the dryer + - ``AgitatorPower``: power to the agitator + - ``AgitatorTorque``: torque resistance for the agitator + - ``AgitatorSpeed``: agitator speed + - ``JacketTemperatureSP``: set point for the jacket heating medium + - ``JacketTemperature``: temperature of the jacket heating medium + - ``DryerTemperatureSP``: set point for the temperature inside the dryer + - ``DryerTemp``: temperature inside the dryer + - ``ClockTime``: sample counter (samples assumed evenly spaced) + + The batches have varying durations, so this dataset is a realistic + candidate for alignment (see :func:`process_improve.batch.batch_dtw` and + :func:`process_improve.batch.resample_to_reference`). + + Returns + ------- + dict[Hashable, pd.DataFrame] + Standard batch-data dictionary: 71 batches, each a dataframe of + 11 numeric columns (10 tags plus ``ClockTime``). + + Source + ------ + Garcia-Munoz, S., "Batch process improvement using latent variable + methods", PhD thesis, McMaster University, 2004. Also analyzed in Wold, + Kettaneh-Wold, MacGregor and Dunn, "Batch Process Modeling and MSPC", + Comprehensive Chemometrics, Elsevier, 2009. + + Examples + -------- + >>> from process_improve.batch.datasets import load_dryer + >>> batches = load_dryer() + >>> "DryerTemp" in next(iter(batches.values())).columns + True + """ + return _load_melted_csv("dryer.csv", batch_id_col="batch_id") + + +def load_batch_fake_data() -> dict[Hashable, pd.DataFrame]: + """Return a small synthetic batch dataset. + + Simulated trajectory data for quick examples and tests: two temperature + tags and one pressure tag per batch, plus ``UCI_minutes`` (minutes since + the start of the batch). The wall-clock timestamp column in the raw CSV + is dropped, so all returned columns are numeric. + + Returns + ------- + dict[Hashable, pd.DataFrame] + Standard batch-data dictionary of synthetic batches, each a dataframe + with columns ``UCI_minutes``, ``Temp1``, ``Temp2``, and ``Pressure1``. + + Examples + -------- + >>> from process_improve.batch.datasets import load_batch_fake_data + >>> batches = load_batch_fake_data() + >>> sorted(next(iter(batches.values())).columns) + ['Pressure1', 'Temp1', 'Temp2', 'UCI_minutes'] + """ + return _load_melted_csv("batch-fake-data.csv", batch_id_col="Batch", drop_columns=["DateTime"]) diff --git a/tests/batch/test_datasets.py b/tests/batch/test_datasets.py new file mode 100644 index 00000000..48de21e0 --- /dev/null +++ b/tests/batch/test_datasets.py @@ -0,0 +1,45 @@ +"""Tests for the bundled batch dataset loaders.""" + +import pandas as pd + +from process_improve.batch.data_input import check_valid_batch_dict +from process_improve.batch.datasets import load_batch_fake_data, load_dryer, load_nylon + + +def test_load_nylon() -> None: + """Nylon loads as a valid batch dict: 57 batches of 10 numeric tags, no NaN.""" + batches = load_nylon() + assert len(batches) == 57 + assert check_valid_batch_dict(batches, no_nan=True) + first = next(iter(batches.values())) + assert list(first.columns) == [f"Tag{i:02d}" for i in range(1, 11)] + assert "batch_id" not in first.columns + + +def test_load_dryer() -> None: + """Dryer loads as a valid batch dict with the documented tag names.""" + batches = load_dryer() + assert len(batches) == 71 + assert check_valid_batch_dict(batches) + first = next(iter(batches.values())) + assert "DryerTemp" in first.columns + assert "ClockTime" in first.columns + assert "batch_id" not in first.columns + + +def test_load_batch_fake_data() -> None: + """Synthetic data loads with only numeric columns (DateTime dropped).""" + batches = load_batch_fake_data() + assert len(batches) >= 2 + assert check_valid_batch_dict(batches) + first = next(iter(batches.values())) + assert "DateTime" not in first.columns + assert "Batch" not in first.columns + + +def test_loaders_reset_row_index() -> None: + """Each per-batch frame counts its own samples from zero.""" + for batches in (load_nylon(), load_dryer()): + for batch in batches.values(): + assert isinstance(batch.index, pd.RangeIndex) + assert batch.index[0] == 0 From 1b64df60a26f0b1a04db862a233592054ee602d9 Mon Sep 17 00:00:00 2001 From: Claude Date: Thu, 16 Jul 2026 21:07:15 +0000 Subject: [PATCH 3/8] Complete the melted/wide/dict converter stubs in batch data_input melted_to_wide returned an empty dict, wide_to_melted an empty frame, and wide_to_dict None; the tests pinned that placeholder behaviour. All three are now implemented by composing the existing, tested conversions (melted_to_dict, dict_to_wide, dict_to_melted), and wide_to_dict inverts the (tag, sequence) pivot directly. The tests now assert round-trip identities instead of emptiness, including a truncation-aligned round-trip on the real nylon dataset. wide_to_dict and the new dataset loaders are exported from the batch namespace. Co-Authored-By: Claude Fable 5 Claude-Session: https://claude.ai/code/session_01SyEpexCyqHn1rwTMYQSiP6 --- src/process_improve/batch/__init__.py | 11 +++ src/process_improve/batch/data_input.py | 100 +++++++++++++++++++----- tests/batch/test_data_input.py | 74 +++++++++++++++--- 3 files changed, 153 insertions(+), 32 deletions(-) diff --git a/src/process_improve/batch/__init__.py b/src/process_improve/batch/__init__.py index 5bb3e272..f7bfd780 100644 --- a/src/process_improve/batch/__init__.py +++ b/src/process_improve/batch/__init__.py @@ -6,8 +6,14 @@ dict_to_wide, melted_to_dict, melted_to_wide, + wide_to_dict, wide_to_melted, ) +from process_improve.batch.datasets import ( + load_batch_fake_data, + load_dryer, + load_nylon, +) from process_improve.batch.features import ( cross, f_agemax, @@ -63,8 +69,13 @@ "f_std", "f_sum", "find_reference_batch", + # Dataset loaders + "load_batch_fake_data", + "load_dryer", + "load_nylon", "melted_to_dict", "melted_to_wide", "resample_to_reference", + "wide_to_dict", "wide_to_melted", ] diff --git a/src/process_improve/batch/data_input.py b/src/process_improve/batch/data_input.py index 943f9e4e..d1529e33 100644 --- a/src/process_improve/batch/data_input.py +++ b/src/process_improve/batch/data_input.py @@ -170,34 +170,92 @@ def melted_to_dict(in_df: pd.DataFrame, batch_id_col: str) -> dict: return {batch_id: batch for batch_id, batch in in_df.groupby(batch_id_col)} # noqa: C416 -def melted_to_wide(in_df: pd.DataFrame, batch_id_col: str) -> dict: - """Convert aligned melted data to wide format.""" +def melted_to_wide(in_df: pd.DataFrame, batch_id_col: str, group_by_batch: bool = False) -> pd.DataFrame: + """ + Convert aligned melted data to wide format. + + Parameters + ---------- + in_df : pd.DataFrame + Melted batch data: all batches stacked vertically, one column per tag, + with a batch-identifier column. The batches must be aligned (the same + number of rows per batch), because the wide format is only meaningful + for aligned data. + batch_id_col : str + Name of the column holding the batch identifier. + group_by_batch : bool, optional + Passed through to :func:`dict_to_wide`; controls whether the 2-level + column index is ordered ``(tag, sequence)`` (default) or + ``(sequence, tag)``. + + Returns + ------- + pd.DataFrame + Wide-format dataframe: one row per batch, 2-level column index. + """ if batch_id_col not in in_df: raise ValueError(f"The `batch_id_col` column {batch_id_col!r} does not exist in the incoming dataframe.") - return {} - - # max_places = int(np.ceil(np.log10(aligned_df["_sequence_"].max()))) - # aligned_wide_df = aligned_df.pivot(index="batch_id", columns="_sequence_") - # new_labels = [ - # "-".join(item) - # for item in zip( - # aligned_wide_df.columns.get_level_values(0), - # [str(val).zfill(max_places) for val in aligned_wide_df.columns.get_level_values(1)], - # ) - # ] - # aligned_wide_df.columns = new_labels - # TODO: add the column multilevel column index. - # return dict_to_wide(melted_to_dict(in_df, batch_id_col)) + return dict_to_wide(melted_to_dict(in_df, batch_id_col), group_by_batch=group_by_batch) + + +def wide_to_dict(in_df: pd.DataFrame) -> dict: + """ + Convert wide-format batch data back to the standard dict format. + + Inverts :func:`dict_to_wide`: each row of the wide frame becomes one entry + in the dictionary, with the 2-level ``(tag, sequence)`` column index + pivoted back to a per-batch dataframe of one column per tag, indexed by + sequence. Accepts either column-level ordering (``(tag, sequence)`` or the + ``group_by_batch=True`` variant ``(sequence, tag)``). + + Parameters + ---------- + in_df : pd.DataFrame + Wide-format batch data: one row per batch, 2-level column index with + levels named ``tag`` and ``sequence``. + + Returns + ------- + dict + Standard batch-data dictionary: keys are the wide frame's row index + (batch identifiers), values are per-batch dataframes. + """ + if in_df.columns.nlevels != 2 or set(in_df.columns.names) != {"tag", "sequence"}: + raise ValueError( + "The wide dataframe must have a 2-level column index with levels named " + f"'tag' and 'sequence'; got levels {in_df.columns.names}." + ) + out = {} + for batch_id, row in in_df.iterrows(): + # Series.unstack is the direct inverse of the pivot in `dict_to_wide`; + # pivot_table would aggregate and lose the (tag, sequence) fidelity. + batch = row.unstack(level="tag") # noqa: PD010 + batch.index.name = None + batch.columns.name = None + out[batch_id] = batch + return out def wide_to_melted(in_df: pd.DataFrame) -> pd.DataFrame: - """Convert wide-format batch data to melted format. Not yet implemented.""" - # dict_to_melted(dict_to_wide(in_df)) - return pd.DataFrame() + """ + Convert wide-format batch data to melted format. + Inverts the melted-to-wide direction: the wide frame (one row per batch, + 2-level column index) is expanded back to a melted frame with all batches + stacked vertically, one column per tag, and a ``batch_id`` column. -def wide_to_dict() -> None: - """Convert wide-format batch data to dict format. Not yet implemented.""" + Parameters + ---------- + in_df : pd.DataFrame + Wide-format batch data: one row per batch, 2-level column index with + levels named ``tag`` and ``sequence``. + + Returns + ------- + pd.DataFrame + Melted batch data with a ``batch_id`` column. + """ + return dict_to_melted(wide_to_dict(in_df), insert_batch_id_column=True) def melt_df_to_series(in_df: pd.DataFrame, exclude_columns: list | None = None, name: str | None = None) -> pd.Series: diff --git a/tests/batch/test_data_input.py b/tests/batch/test_data_input.py index 97ab18e4..06130245 100644 --- a/tests/batch/test_data_input.py +++ b/tests/batch/test_data_input.py @@ -9,6 +9,7 @@ melt_df_to_series, melted_to_dict, melted_to_wide, + wide_to_dict, wide_to_melted, ) @@ -36,21 +37,72 @@ def test_melted_to_dict(nylon_raw_melteddata: pd.DataFrame) -> None: assert len(out) == 57 -def test_melted_to_wide(nylon_raw_melteddata: pd.DataFrame) -> None: - """Test conversion from melted format to wide format.""" - _ = melted_to_wide(nylon_raw_melteddata, batch_id_col="batch_id") - # assert out.shape == pytest.approx([2, 3]) +def test_melted_to_wide(aligned_batch_dict: dict) -> None: + """melted_to_wide matches the dict_to_wide(melted_to_dict(...)) composition.""" + melted = dict_to_melted(aligned_batch_dict) + out = melted_to_wide(melted, batch_id_col="batch_id") + expected = dict_to_wide(aligned_batch_dict) + pd.testing.assert_frame_equal(out, expected) + assert list(out.columns.names) == ["tag", "sequence"] + + +def test_melted_to_wide_missing_batch_id_column() -> None: + """A missing batch-identifier column raises a clear error.""" + with pytest.raises(ValueError, match="does not exist"): + melted_to_wide(pd.DataFrame({"x": [1, 2]}), batch_id_col="batch_id") + +def test_wide_to_dict_round_trip(aligned_batch_dict: dict) -> None: + """wide_to_dict inverts dict_to_wide.""" + out = wide_to_dict(dict_to_wide(aligned_batch_dict)) + assert set(out.keys()) == set(aligned_batch_dict.keys()) + for batch_id, original in aligned_batch_dict.items(): + recovered = out[batch_id] + pd.testing.assert_frame_equal( + recovered.sort_index(axis=1), original.sort_index(axis=1), check_index_type=False + ) -def test_wide_to_melted() -> None: - """Test conversion from wide format to melted format.""" - out = wide_to_melted(pd.DataFrame({"x": [1, 2]})) - assert isinstance(out, pd.DataFrame) - assert out.empty +def test_wide_to_dict_group_by_batch_round_trip(aligned_batch_dict: dict) -> None: + """wide_to_dict also inverts the group_by_batch=True column ordering.""" + out = wide_to_dict(dict_to_wide(aligned_batch_dict, group_by_batch=True)) + for batch_id, original in aligned_batch_dict.items(): + pd.testing.assert_frame_equal( + out[batch_id].sort_index(axis=1), original.sort_index(axis=1), check_index_type=False + ) -def test_wide_to_dict() -> None: - """Test conversion from wide format to dictionary.""" + +def test_wide_to_dict_rejects_flat_columns() -> None: + """A frame without the 2-level (tag, sequence) column index is rejected.""" + with pytest.raises(ValueError, match="2-level column index"): + wide_to_dict(pd.DataFrame({"x": [1, 2]})) + + +def test_wide_to_melted_round_trip(aligned_batch_dict: dict) -> None: + """wide_to_melted recovers a melted frame equivalent to dict_to_melted's.""" + out = wide_to_melted(dict_to_wide(aligned_batch_dict)) + expected = dict_to_melted(aligned_batch_dict) + assert "batch_id" in out.columns + assert out.shape[0] == expected.shape[0] + for batch_id in aligned_batch_dict: + recovered = out[out["batch_id"] == batch_id].drop(columns="batch_id").reset_index(drop=True) + original = expected[expected["batch_id"] == batch_id].drop(columns="batch_id").reset_index(drop=True) + pd.testing.assert_frame_equal( + recovered.sort_index(axis=1), original.sort_index(axis=1), check_index_type=False + ) + + +def test_melted_to_wide_on_real_data(nylon_raw_melteddata: pd.DataFrame) -> None: + """Wide conversion of (truncation-aligned) nylon data has one row per batch.""" + # The bundled nylon batches differ slightly in duration; truncate to the + # shortest so the wide format is meaningful. + batches = melted_to_dict(nylon_raw_melteddata, batch_id_col="batch_id") + shortest = min(batch.shape[0] for batch in batches.values()) + truncated = {k: v.iloc[:shortest].reset_index(drop=True) for k, v in batches.items()} + wide = dict_to_wide({k: v.drop(columns="batch_id") for k, v in truncated.items()}) + assert wide.shape == (57, 10 * shortest) + recovered = wide_to_dict(wide) + assert len(recovered) == 57 def test_dict_to_melted_default_inserts_batch_id(aligned_batch_dict: dict) -> None: From 76dcdda30ba4137c32258431ee19a808be786002 Mon Sep 17 00:00:00 2001 From: Claude Date: Thu, 16 Jul 2026 21:09:31 +0000 Subject: [PATCH 4/8] Sweep lint and dead-code debt in the batch subpackage - plotting.py: fix the silently-ignored 'resonsive' Plotly config typo; replace stdlib random.seed/shuffle with the package-wide check_random_state contract; replace the inert file-level 'flake8: noqa: C901' directive with per-function noqa codes (it was blanket-suppressing C901 for the whole file under ruff). - preprocessing.py: remove the dead guarded plotly import (go was never used in this module) and route the DTW progress print through the module logger. Co-Authored-By: Claude Fable 5 Claude-Session: https://claude.ai/code/session_01SyEpexCyqHn1rwTMYQSiP6 --- src/process_improve/batch/plotting.py | 18 +++++++++--------- src/process_improve/batch/preprocessing.py | 8 +------- 2 files changed, 10 insertions(+), 16 deletions(-) diff --git a/src/process_improve/batch/plotting.py b/src/process_improve/batch/plotting.py index 68e30712..fa3fe512 100644 --- a/src/process_improve/batch/plotting.py +++ b/src/process_improve/batch/plotting.py @@ -3,7 +3,6 @@ # Built-in libraries import json import math -import random from collections.abc import Callable, Sequence from functools import partial from typing import Any, cast @@ -24,6 +23,7 @@ sns = _MissingExtra("seaborn", "plotting") # type: ignore[assignment] plotoffline = _MissingExtra("plotly", "plotting") # type: ignore[assignment] +from .._random import check_random_state from .data_input import check_valid_batch_dict @@ -88,7 +88,7 @@ def plot_to_HTML(filename: str, fig: dict) -> str: # noqa: N802 editable=False, displaylogo=False, showLink=False, - resonsive=True, + responsive=True, ) return plotoffline( figure_or_data=fig, @@ -100,7 +100,7 @@ def plot_to_HTML(filename: str, fig: dict) -> str: # noqa: N802 ) -def plot_all_batches_per_tag( # noqa: PLR0912, PLR0913 +def plot_all_batches_per_tag( # noqa: C901, PLR0912, PLR0913 df_dict: dict, tag: str, tag_y2: str | None = None, @@ -171,9 +171,10 @@ def plot_all_batches_per_tag( # noqa: PLR0912, PLR0913 default_line_width = 2 unique_items = list(df_dict.keys()) n_colours = len(unique_items) - random.seed(13) colours = list(sns.husl_palette(n_colours)) - random.shuffle(colours) + # Deterministic shuffle so adjacent batches get visually distinct colours, + # reproducibly across calls. Uses the package-wide RNG contract. + check_random_state(13).shuffle(colours) colours = [get_rgba_from_triplet(c, as_string=True) for c in colours] line_styles = {k: dict(width=default_line_width, color=v) for k, v in zip(unique_items, colours, strict=False)} for key, val in batches_to_highlight.items(): @@ -305,10 +306,10 @@ def colours_per_batch_id( """ if colour_map is None: colour_map = partial(sns.color_palette, "hls") - random.seed(13) n_colours = len(batch_ids) colours = list(colour_map(n_colours)) if not (use_default_colour) else [(0.5, 0.5, 0.5)] * n_colours - random.shuffle(colours) + # Deterministic shuffle (package-wide RNG contract) for reproducible colours. + check_random_state(13).shuffle(colours) colours = [get_rgba_from_triplet(c, as_string=True) for c in colours] colour_assignment = { key: dict(width=default_line_width, color=val) @@ -329,8 +330,7 @@ def colours_per_batch_id( return colour_assignment -# flake8: noqa: C901 -def plot_multitags( # noqa: PLR0912, PLR0913, PLR0915 +def plot_multitags( # noqa: C901, PLR0912, PLR0913, PLR0915 df_dict: dict, batch_list: list | None = None, tag_list: list | None = None, diff --git a/src/process_improve/batch/preprocessing.py b/src/process_improve/batch/preprocessing.py index fbbc84e8..735514ed 100644 --- a/src/process_improve/batch/preprocessing.py +++ b/src/process_improve/batch/preprocessing.py @@ -8,12 +8,6 @@ import scipy as sp from tqdm import tqdm -try: - import plotly.graph_objects as go -except ImportError: # pragma: no cover - exercised via env-without-plotly - from process_improve._extras import _MissingExtra - go = _MissingExtra("plotly", "plotting") # type: ignore[assignment] - from ..multivariate.methods import PCA, MCUVScaler from .alignment_helpers import backtrack_optimal_path, distance_matrix from .data_input import check_valid_batch_dict, dict_to_wide, melted_to_dict @@ -366,7 +360,7 @@ def batch_dtw( # noqa: C901, PLR0915 iter_step = 0 while (np.linalg.norm(delta_weight) > settings["tolerance"]) and (iter_step <= settings["maximum_iterations"]): if settings["show_progress"]: - print(f"Iter = {iter_step} and norm = {np.linalg.norm(delta_weight)}") # noqa: T201 + logger.info("Iter = %d and norm = %g", iter_step, float(np.linalg.norm(delta_weight))) iter_step += 1 logger.debug( From 68479e03412d971421570f3bb9c9de17cefaa87f Mon Sep 17 00:00:00 2001 From: Claude Date: Thu, 16 Jul 2026 21:15:47 +0000 Subject: [PATCH 5/8] Add BatchPCA: batchwise-unfolded (multiway) PCA for batch data BatchPCA brings the Nomikos-MacGregor batch modelling approach to the library. It unfolds an aligned batch-data dictionary batchwise (one row per batch) via dict_to_wide, optionally joins an initial-conditions (Z) block onto that row, mean-centres and scales every column with MCUVScaler (which removes the average trajectory), and fits the existing multivariate.PCA by composition. Batch-level scores, SPE and Hotelling's T2 with control limits, plus contribution and score/SPE/loading plot forwarders, are exposed as batch-indexed views. Loadings keep the (tag, sequence) index so they reshape to a variable-by-time grid. Until now the batch module stopped at feature extraction with no multiway model; the unfold-then-PCA path existed only privately inside find_reference_batch. This is the headline capability the batch case study in the book depends on. Co-Authored-By: Claude Fable 5 Claude-Session: https://claude.ai/code/session_01SyEpexCyqHn1rwTMYQSiP6 --- src/process_improve/batch/__init__.py | 5 +- src/process_improve/batch/_batch_pca.py | 404 ++++++++++++++++++++++++ tests/batch/test_batch_pca.py | 127 ++++++++ 3 files changed, 535 insertions(+), 1 deletion(-) create mode 100644 src/process_improve/batch/_batch_pca.py create mode 100644 tests/batch/test_batch_pca.py diff --git a/src/process_improve/batch/__init__.py b/src/process_improve/batch/__init__.py index f7bfd780..3d17d207 100644 --- a/src/process_improve/batch/__init__.py +++ b/src/process_improve/batch/__init__.py @@ -1,5 +1,6 @@ -"""Batch process data analysis: alignment, feature extraction, and visualization.""" +"""Batch process data analysis: alignment, feature extraction, modelling, and visualization.""" +from process_improve.batch._batch_pca import BatchPCA from process_improve.batch.data_input import ( check_valid_batch_dict, dict_to_melted, @@ -42,6 +43,8 @@ ) __all__ = [ + # Modelling + "BatchPCA", # Preprocessing/alignment "batch_dtw", # Data input/conversion diff --git a/src/process_improve/batch/_batch_pca.py b/src/process_improve/batch/_batch_pca.py new file mode 100644 index 00000000..2be9dd43 --- /dev/null +++ b/src/process_improve/batch/_batch_pca.py @@ -0,0 +1,404 @@ +# (c) Kevin Dunn, 2010-2026. MIT License. Based on own private work over the years. +"""Batchwise-unfolded (multiway) PCA for batch trajectory data. + +Implements the batch modelling approach of Nomikos and MacGregor: align the +batches, unfold them batchwise (one row per batch), mean-centre and scale each +column (which removes the average trajectory), and fit an ordinary PCA on the +result. The fitted model summarizes the deviations of every batch from the +average trajectory, so batches can be compared, and new batches diagnosed, +in a low-dimensional score space with Hotelling's T2 and SPE limits. + +Initial conditions (the Z block: one row of pre-batch measurements per batch) +can be appended to the unfolded row, so the model sees ``[Z | X-unfolded]`` +with exactly one row per batch. + +See Wold, Kettaneh-Wold, MacGregor and Dunn, "Batch Process Modeling and +MSPC", Comprehensive Chemometrics, Elsevier, 2009, for the methodology. +""" + +from __future__ import annotations + +import functools +import typing + +import pandas as pd +from sklearn.base import BaseEstimator, TransformerMixin +from sklearn.utils.validation import check_is_fitted + +from ..multivariate._diagnostics import spe_contributions as _spe_contributions +from ..multivariate._diagnostics import t2_contributions as _t2_contributions +from ..multivariate._limits import score_limit as _score_limit +from ..multivariate._limits import spe_limit as _spe_limit +from ..multivariate._pca import PCA +from ..multivariate._preprocessing import MCUVScaler +from ..multivariate.plots import explained_variance_plot as _explained_variance_plot +from ..multivariate.plots import loading_plot as _loading_plot +from ..multivariate.plots import score_plot as _score_plot +from ..multivariate.plots import spe_plot as _spe_plot +from ..multivariate.plots import t2_plot as _t2_plot +from .data_input import check_valid_batch_dict, dict_to_wide + +if typing.TYPE_CHECKING: + from collections.abc import Callable, Hashable + + import numpy as np + from sklearn.utils import Bunch + + +def _pca_method(fn: Callable[..., typing.Any]) -> Callable[..., typing.Any]: + """Wrap a module-level ``fn(model, ...)`` as a method forwarding ``self._pca``. + + The convenience plots, limits, and contribution functions in + :mod:`process_improve.multivariate` read only fitted attributes that the + internal PCA model carries (``scores_``, ``loadings_``, ``spe_``, + ``scaling_factor_for_scores_``, ``n_components``, ``n_samples_``), so the + :class:`BatchPCA` methods forward to them with the wrapped estimator as + the ``model`` argument. Mirrors ``_model_method`` (ENG-05), so ``help`` + and ``inspect.signature`` report the underlying function. + """ + + @functools.wraps(fn) + def method(self: BatchPCA, *args, **kwargs) -> object: + check_is_fitted(self, "loadings_") + return fn(self._pca, *args, **kwargs) + + return method + + +class BatchPCA(TransformerMixin, BaseEstimator): + """Batchwise-unfolded (multiway) PCA on aligned batch trajectory data. + + Each batch becomes one row of the model matrix: the trajectories are + unfolded batchwise via :func:`process_improve.batch.dict_to_wide`, the + optional initial-conditions block is joined on, every column is + mean-centred and (optionally) scaled to unit variance with + :class:`process_improve.multivariate.MCUVScaler`, and an ordinary + :class:`process_improve.multivariate.PCA` is fitted to the result. + Centring the unfolded columns removes the average trajectory, so the + components model the batch-to-batch deviations. + + The batches must be aligned before fitting: every batch must have the + same number of samples (see + :func:`process_improve.batch.resample_to_reference` and + :func:`process_improve.batch.batch_dtw`), and no missing values are + allowed. + + Parameters + ---------- + n_components : int + Number of principal components to extract. + scale : bool, default=True + Scale each unfolded column to unit variance after centring. Centring + always happens (it removes the average trajectory); set this to False + to keep the columns in their centred, unscaled units. + group_by_batch : bool, default=False + Ordering of the unfolded column index, passed to + :func:`process_improve.batch.dict_to_wide`: ``False`` groups all time + samples of a tag together (``(tag, sequence)``); ``True`` groups all + tags of a time sample together (``(sequence, tag)``). + algorithm : str, default="auto" + Fitting algorithm, passed to + :class:`process_improve.multivariate.PCA`. + + Attributes (after fitting) + -------------------------- + scores_ : pd.DataFrame of shape (n_batches, n_components) + Batch-level scores; one row per batch, indexed by batch identifier. + loadings_ : pd.DataFrame of shape (n_unfolded_features, n_components) + Loadings, indexed by the 2-level unfolded column index, so the + trajectory part reshapes to a (tag, time) grid. Initial-condition + rows (if any) carry an empty string in the ``sequence`` level. + spe_ : pd.DataFrame of shape (n_batches, n_components) + Per-batch SPE after each component (residual scale). + hotellings_t2_ : pd.DataFrame of shape (n_batches, n_components) + Per-batch cumulative Hotelling's T2. + explained_variance_ : np.ndarray of shape (n_components,) + Variance explained by each component. + r2_per_component_, r2_cumulative_ : pd.Series of length n_components + Fractional and cumulative R2 of the unfolded matrix. + n_batches_ : int + Number of batches in the training set. + n_tags_ : int + Number of trajectory tags per batch. + n_timesteps_ : int + Number of (aligned) time samples per batch. + n_initial_conditions_ : int + Number of initial-condition (Z) columns; zero when none were given. + batch_ids_ : list + Batch identifiers, in model-row order. + tag_names_ : list + Trajectory tag names. + initial_condition_names_ : list + Initial-condition column names (empty when none were given). + time_index_ : list + The aligned sequence values (0, 1, ..., ``n_timesteps_`` - 1). + center_, scale_ : pd.Series of length n_unfolded_features + The per-column centring and scaling applied before the PCA fit. + + Examples + -------- + >>> from process_improve.batch import BatchPCA, load_nylon, resample_to_reference + >>> batches = load_nylon() + >>> tags = list(next(iter(batches.values())).columns) + >>> aligned = resample_to_reference(batches, columns_to_align=tags, reference_batch="1") + >>> model = BatchPCA(n_components=3).fit(aligned) + >>> model.scores_.shape + (57, 3) + + See Also + -------- + process_improve.multivariate.PCA : the underlying estimator. + + References + ---------- + Nomikos, P. and MacGregor, J.F., "Monitoring of Batch Processes Using + Multi-Way Principal Component Analysis", AIChE Journal, 40, 1361-1375, + 1994. + + Wold, S., Kettaneh-Wold, N., MacGregor, J.F. and Dunn, K.G., "Batch + Process Modeling and MSPC", Comprehensive Chemometrics, Elsevier, 2009. + """ + + _parameter_constraints: typing.ClassVar = { + "n_components": [int, None], + "scale": [bool], + "group_by_batch": [bool], + "algorithm": [str], + } + + def __init__( + self, + n_components: int, + *, + scale: bool = True, + group_by_batch: bool = False, + algorithm: str = "auto", + ) -> None: + self.n_components = n_components + self.scale = scale + self.group_by_batch = group_by_batch + self.algorithm = algorithm + + def _unfold( + self, + batches: dict[Hashable, pd.DataFrame], + initial_conditions: pd.DataFrame | None, + ) -> pd.DataFrame: + """Unfold the batches batchwise and join the initial-conditions block. + + Returns the one-row-per-batch ``[Z | X-unfolded]`` matrix with a + 2-level column index throughout: trajectory columns keep their + ``(tag, sequence)`` labels; initial-condition columns are labelled + ``(name, "")`` since they carry no time axis. + """ + check_valid_batch_dict(batches, no_nan=True) + wide = dict_to_wide(batches, group_by_batch=self.group_by_batch) + if initial_conditions is None: + return wide + + if not isinstance(initial_conditions, pd.DataFrame): + raise TypeError( + "initial_conditions must be a pandas DataFrame indexed by batch identifier; " + f"got {type(initial_conditions).__name__}." + ) + if set(initial_conditions.index) != set(wide.index): + missing = set(wide.index) - set(initial_conditions.index) + extra = set(initial_conditions.index) - set(wide.index) + raise ValueError( + "initial_conditions must have exactly one row per batch. " + f"Missing batch ids: {sorted(missing, key=str)}; unmatched extra ids: {sorted(extra, key=str)}." + ) + z_wide = initial_conditions.reindex(wide.index) + if z_wide.select_dtypes(include="number").shape[1] != z_wide.shape[1]: + raise ValueError("All initial_conditions columns must be numeric.") + if z_wide.isna().to_numpy().sum() > 0: + raise ValueError("No missing values allowed in initial_conditions.") + if self.group_by_batch: + tuples = [("", name) for name in z_wide.columns] + else: + tuples = [(name, "") for name in z_wide.columns] + z_wide.columns = pd.MultiIndex.from_tuples(tuples, names=wide.columns.names) + return pd.concat([z_wide, wide], axis=1) + + def fit( + self, + X: dict[Hashable, pd.DataFrame], + y: object = None, # noqa: ARG002 + *, + initial_conditions: pd.DataFrame | None = None, + ) -> BatchPCA: + """Fit the batchwise-unfolded PCA model. + + Parameters + ---------- + X : dict[Hashable, pd.DataFrame] + Standard batch-data dictionary of aligned batches: keys are batch + identifiers, values are per-batch dataframes with identical + all-numeric columns and the same number of rows. No missing + values. + y : ignored + Present for sklearn Pipeline compatibility. + initial_conditions : pd.DataFrame, optional + The Z block: one row per batch (indexed by the same batch + identifiers as ``X``), one column per pre-batch measurement. + Joined onto the unfolded row before centring and scaling. + + Returns + ------- + self : BatchPCA + """ + wide = self._unfold(X, initial_conditions) + + scaler = MCUVScaler().fit(wide) + if not self.scale: + scaler.scale_ = pd.Series(1.0, index=scaler.scale_.index) + mcuv = pd.DataFrame(scaler.transform(wide).to_numpy(), index=wide.index, columns=wide.columns) + + self._scaler = scaler + self._pca = PCA(n_components=self.n_components, algorithm=self.algorithm).fit(mcuv) + + # Batch-shaped views over the fitted internal model. The internal PCA + # was fitted on the wide frame directly, so its row index is already + # the batch identifiers and its feature index is the 2-level unfolded + # column index. + self.feature_columns_ = wide.columns + self.scores_ = self._pca.scores_ + self.loadings_ = self._pca.loadings_ + self.spe_ = self._pca.spe_ + self.hotellings_t2_ = self._pca.hotellings_t2_ + self.explained_variance_ = self._pca.explained_variance_ + self.r2_per_component_ = self._pca.r2_per_component_ + self.r2_cumulative_ = self._pca.r2_cumulative_ + self.scaling_factor_for_scores_ = self._pca.scaling_factor_for_scores_ + + first_batch = X[next(iter(X.keys()))] + self.batch_ids_ = list(wide.index) + self.n_batches_ = len(self.batch_ids_) + self.tag_names_ = list(first_batch.columns) + self.n_tags_ = len(self.tag_names_) + self.n_timesteps_ = int(first_batch.shape[0]) + self.time_index_ = list(range(self.n_timesteps_)) + if initial_conditions is None: + self.initial_condition_names_ = [] + self.n_initial_conditions_ = 0 + else: + self.initial_condition_names_ = list(initial_conditions.columns) + self.n_initial_conditions_ = len(self.initial_condition_names_) + self.center_ = scaler.center_ + self.scale_ = scaler.scale_ + self.n_samples_ = self._pca.n_samples_ + return self + + def _scaled_wide( + self, + batches: dict[Hashable, pd.DataFrame], + initial_conditions: pd.DataFrame | None, + ) -> pd.DataFrame: + """Unfold new batches and apply the training centring/scaling.""" + check_is_fitted(self, "loadings_") + wide = self._unfold(batches, initial_conditions) + if list(wide.columns) != list(self.feature_columns_): + raise ValueError( + "The new batches do not unfold to the training column layout. " + f"Expected {len(self.feature_columns_)} unfolded columns " + f"({self.n_tags_} tags x {self.n_timesteps_} samples" + + (f" + {self.n_initial_conditions_} initial conditions" if self.n_initial_conditions_ else "") + + f"); got {len(wide.columns)}. Align new batches to the training " + "length and pass the same tags and initial-condition columns." + ) + return pd.DataFrame(self._scaler.transform(wide).to_numpy(), index=wide.index, columns=wide.columns) + + def transform( + self, + X: dict[Hashable, pd.DataFrame], + *, + initial_conditions: pd.DataFrame | None = None, + ) -> pd.DataFrame: + """Project new (complete, aligned) batches onto the model. + + Parameters + ---------- + X : dict[Hashable, pd.DataFrame] + Standard batch-data dictionary of aligned batches with the same + tags and number of samples as the training data. + initial_conditions : pd.DataFrame, optional + The Z block for the new batches; required if (and only if) the + model was fitted with one. + + Returns + ------- + pd.DataFrame of shape (n_new_batches, n_components) + Batch-level scores, indexed by batch identifier. + """ + return self._pca.transform(self._scaled_wide(X, initial_conditions)) + + def fit_transform( + self, + X: dict[Hashable, pd.DataFrame], + y: object = None, + *, + initial_conditions: pd.DataFrame | None = None, + ) -> pd.DataFrame: + """Fit the model and return the training batch scores.""" + self.fit(X, y, initial_conditions=initial_conditions) + return self.scores_ + + def diagnose( + self, + X: dict[Hashable, pd.DataFrame], + *, + initial_conditions: pd.DataFrame | None = None, + ) -> Bunch: + """Project new batches and compute their monitoring diagnostics. + + Parameters + ---------- + X : dict[Hashable, pd.DataFrame] + Standard batch-data dictionary of aligned batches with the same + tags and number of samples as the training data. + initial_conditions : pd.DataFrame, optional + The Z block for the new batches; required if (and only if) the + model was fitted with one. + + Returns + ------- + result : sklearn.utils.Bunch + With keys ``scores`` (DataFrame, one row per batch), + ``hotellings_t2`` (DataFrame, cumulative per component), and + ``spe`` (Series). Compare against :meth:`hotellings_t2_limit` + and :meth:`spe_limit` to flag abnormal batches. + """ + return self._pca.diagnose(self._scaled_wide(X, initial_conditions)) + + def hotellings_t2_limit(self, conf_level: float = 0.95) -> float: + """Hotelling's T2 limit at the given confidence level.""" + check_is_fitted(self, "loadings_") + return self._pca.hotellings_t2_limit(conf_level=conf_level) + + def ellipse_coordinates( + self, + score_horiz: int, + score_vert: int, + conf_level: float = 0.95, + n_points: int = 100, + ) -> tuple[np.ndarray, np.ndarray]: + """Coordinates of the T2 confidence ellipse for a score plot.""" + check_is_fitted(self, "loadings_") + return self._pca.ellipse_coordinates( + score_horiz=score_horiz, + score_vert=score_vert, + conf_level=conf_level, + n_points=n_points, + ) + + # Convenience methods forwarding to the standalone multivariate functions + # with the internal (batchwise-unfolded) PCA as the model argument. + score_plot = _pca_method(_score_plot) + spe_plot = _pca_method(_spe_plot) + t2_plot = _pca_method(_t2_plot) + loading_plot = _pca_method(_loading_plot) + explained_variance_plot = _pca_method(_explained_variance_plot) + spe_limit = _pca_method(_spe_limit) + score_limit = _pca_method(_score_limit) + t2_contributions = _pca_method(_t2_contributions) + spe_contributions = _pca_method(_spe_contributions) diff --git a/tests/batch/test_batch_pca.py b/tests/batch/test_batch_pca.py new file mode 100644 index 00000000..2f0290ff --- /dev/null +++ b/tests/batch/test_batch_pca.py @@ -0,0 +1,127 @@ +"""Tests for BatchPCA: batchwise-unfolded (multiway) PCA on batch data.""" + +import numpy as np +import pandas as pd +import pytest + +from process_improve.batch._batch_pca import BatchPCA +from process_improve.batch.datasets import load_nylon +from process_improve.batch.preprocessing import resample_to_reference + + +@pytest.fixture +def aligned_nylon() -> dict: + """Nylon batches resampled to a common length, ready for unfolding.""" + batches = load_nylon() + tags = list(next(iter(batches.values())).columns) + return resample_to_reference(batches, columns_to_align=tags, reference_batch=1) + + +def _synthetic_good_batches(n_batches: int = 12, n_timesteps: int = 40, *, seed_offset: int = 0) -> dict: + """Build a set of similar, smoothly-varying good batches (three tags).""" + t = np.linspace(0, 1, n_timesteps) + out = {} + for i in range(n_batches): + rng = np.random.default_rng(i + seed_offset) + out[f"g{i}"] = pd.DataFrame( + { + "A": np.sin(2 * np.pi * t) + 0.02 * rng.standard_normal(n_timesteps), + "B": t**2 + 0.02 * rng.standard_normal(n_timesteps), + "C": np.exp(t) + 0.02 * rng.standard_normal(n_timesteps), + } + ) + return out + + +def test_fit_on_real_nylon(aligned_nylon: dict) -> None: + """Fitting on aligned nylon yields batch-level scores and reshapeable loadings.""" + model = BatchPCA(n_components=3).fit(aligned_nylon) + assert model.scores_.shape == (57, 3) + assert model.n_batches_ == 57 + assert model.n_tags_ == 10 + # loadings reshape to a (tag, time) grid + grid = model.loadings_.iloc[:, 0].unstack(level="sequence") # noqa: PD010 + assert grid.shape == (model.n_tags_, model.n_timesteps_) + # cumulative R2 is monotone non-decreasing and within [0, 1] + r2 = model.r2_cumulative_.to_numpy() + assert np.all(np.diff(r2) >= -1e-9) + assert 0.0 <= r2[0] <= r2[-1] <= 1.0 + + +def test_transform_round_trip(aligned_nylon: dict) -> None: + """Transforming the training batches reproduces the fitted scores.""" + model = BatchPCA(n_components=3).fit(aligned_nylon) + scores = model.transform(aligned_nylon) + np.testing.assert_allclose(scores.to_numpy(), model.scores_.to_numpy(), atol=1e-8) + + +def test_diagnose_matches_internal_pca(aligned_nylon: dict) -> None: + """Diagnosing the training data matches the fitted diagnostics.""" + model = BatchPCA(n_components=3).fit(aligned_nylon) + result = model.diagnose(aligned_nylon) + np.testing.assert_allclose(result.scores.to_numpy(), model.scores_.to_numpy(), atol=1e-8) + np.testing.assert_allclose( + result.hotellings_t2.to_numpy(), model.hotellings_t2_.to_numpy(), atol=1e-8 + ) + + +def test_fit_returns_self_and_limits(aligned_nylon: dict) -> None: + """Fitting returns self; the limit helpers return positive finite numbers.""" + model = BatchPCA(n_components=3) + assert model.fit(aligned_nylon) is model + assert model.hotellings_t2_limit(0.95) > 0 + assert model.spe_limit(0.95) > 0 + + +def test_injected_fault_is_flagged() -> None: + """A batch with a large deviation trips the SPE and T2 limits.""" + good = _synthetic_good_batches() + model = BatchPCA(n_components=2).fit(good) + + faulty = _synthetic_good_batches(n_batches=1, seed_offset=99) + only = next(iter(faulty.values())) + only["A"] = only["A"] + 1.5 # sustained upset on tag A + result = model.diagnose({"bad": only}) + + assert float(result.spe.iloc[-1]) > model.spe_limit(0.95) + assert float(result.hotellings_t2.iloc[-1, -1]) > model.hotellings_t2_limit(0.95) + + +def test_initial_conditions_join() -> None: + """The Z block is appended to the unfolded row, one row per batch.""" + good = _synthetic_good_batches() + z = pd.DataFrame( + {"charge": [10.0 + i for i in range(12)], "purity": [0.9 - 0.01 * i for i in range(12)]}, + index=list(good.keys()), + ) + model = BatchPCA(n_components=2).fit(good, initial_conditions=z) + # 3 tags x 40 samples + 2 initial conditions + assert model.loadings_.shape[0] == 3 * 40 + 2 + assert model.n_initial_conditions_ == 2 + z_rows = [label for label in model.loadings_.index if label[1] == ""] + assert {label[0] for label in z_rows} == {"charge", "purity"} + + +def test_initial_conditions_index_mismatch_raises() -> None: + """A Z block whose index does not match the batch ids is rejected.""" + good = _synthetic_good_batches(n_batches=4) + z = pd.DataFrame({"charge": [1.0, 2.0, 3.0]}, index=["g0", "g1", "g2"]) # missing g3 + with pytest.raises(ValueError, match="one row per batch"): + BatchPCA(n_components=2).fit(good, initial_conditions=z) + + +def test_transform_wrong_length_raises(aligned_nylon: dict) -> None: + """Projecting batches of a different length gives a clear error.""" + model = BatchPCA(n_components=2).fit(aligned_nylon) + shorter = {k: v.iloc[:50].reset_index(drop=True) for k, v in aligned_nylon.items()} + with pytest.raises(ValueError, match="training column layout"): + model.transform(shorter) + + +def test_unscaled_keeps_centring_only() -> None: + """scale=False leaves the scaling factors at 1.0 (centring still applied).""" + good = _synthetic_good_batches() + model = BatchPCA(n_components=2, scale=False).fit(good) + np.testing.assert_allclose(model.scale_.to_numpy(), 1.0) + # centring is non-trivial + assert np.any(np.abs(model.center_.to_numpy()) > 1e-6) From d8067bc7370530e1e00c7b25e2a77547c057e658 Mon Sep 17 00:00:00 2001 From: Claude Date: Thu, 16 Jul 2026 21:19:47 +0000 Subject: [PATCH 6/8] Add batch model plots: time-varying loadings and contribution-at-time time_varying_loading_plot draws one component's loadings as a function of time (one trace per tag), reading the batchwise-unfolded model's (tag, sequence) loading grid; initial-condition loadings are shown as a marker group before time zero. contribution_at_time_plot renders the per-tag SPE or T2 contribution at a chosen time sample for one batch, to localize an abnormal event to the responsible variable. Both are plotly, guarded behind the plotting extra, and use the package theme. The batch score/SPE/T2 plots reuse the multivariate plot functions on the internal model, so only these two time-resolved plots are new. Co-Authored-By: Claude Fable 5 Claude-Session: https://claude.ai/code/session_01SyEpexCyqHn1rwTMYQSiP6 --- src/process_improve/batch/__init__.py | 6 + src/process_improve/batch/_batch_plots.py | 182 ++++++++++++++++++++++ tests/batch/test_batch_plots.py | 76 +++++++++ 3 files changed, 264 insertions(+) create mode 100644 src/process_improve/batch/_batch_plots.py create mode 100644 tests/batch/test_batch_plots.py diff --git a/src/process_improve/batch/__init__.py b/src/process_improve/batch/__init__.py index 3d17d207..4601862a 100644 --- a/src/process_improve/batch/__init__.py +++ b/src/process_improve/batch/__init__.py @@ -1,6 +1,10 @@ """Batch process data analysis: alignment, feature extraction, modelling, and visualization.""" from process_improve.batch._batch_pca import BatchPCA +from process_improve.batch._batch_plots import ( + contribution_at_time_plot, + time_varying_loading_plot, +) from process_improve.batch.data_input import ( check_valid_batch_dict, dict_to_melted, @@ -49,6 +53,7 @@ "batch_dtw", # Data input/conversion "check_valid_batch_dict", + "contribution_at_time_plot", "cross", "determine_scaling", "dict_to_melted", @@ -79,6 +84,7 @@ "melted_to_dict", "melted_to_wide", "resample_to_reference", + "time_varying_loading_plot", "wide_to_dict", "wide_to_melted", ] diff --git a/src/process_improve/batch/_batch_plots.py b/src/process_improve/batch/_batch_plots.py new file mode 100644 index 00000000..9346c62e --- /dev/null +++ b/src/process_improve/batch/_batch_plots.py @@ -0,0 +1,182 @@ +# (c) Kevin Dunn, 2010-2026. MIT License. Based on own private work over the years. +"""Model-level plots for batchwise-unfolded (multiway) PCA of batch data. + +These complement the batch score / SPE / T2 plots, which are inherited from the +multivariate package (they operate on the internal PCA model). Two plots are +specific to the batch (unfolded) structure: + +- :func:`time_varying_loading_plot`: the loadings of one component drawn as a + function of time, one trace per tag, so the reader sees how each variable + contributes to a component over the batch evolution. +- :func:`contribution_at_time_plot`: for one batch, the per-tag contribution to + SPE or Hotelling's T2 at a chosen time sample, drawn as a bar chart to + diagnose which variable drives an abnormal event. +""" + +from __future__ import annotations + +import typing + +import numpy as np +import pandas as pd + +try: + import plotly.graph_objects as go +except ImportError: # pragma: no cover - exercised via env-without-plotly + from process_improve._extras import _MissingExtra + + go = _MissingExtra("plotly", "plotting") # type: ignore[assignment] + +from ..visualization.themes import DEFAULT_THEME, REFERENCE_LINE_COLOR + +if typing.TYPE_CHECKING: + from ._batch_pca import BatchPCA + + +def _split_loadings(model: BatchPCA, component: int) -> tuple[pd.DataFrame, pd.Series]: + """Split one component's loadings into a (tag x time) grid and a Z series. + + Returns the trajectory loadings reshaped to rows = tags, columns = time + (in fitted order), and the initial-condition loadings as a plain Series + (empty when the model was fitted without a Z block). + """ + loading = model.loadings_.iloc[:, component - 1] + sequence = loading.index.get_level_values("sequence") + is_traj = sequence != "" + traj = loading[is_traj] + # Reshape to tags x time, preserving the fitted tag and time order. + grid = traj.unstack(level="sequence") # noqa: PD010 - direct inverse of the unfold; pivot_table would aggregate + grid = grid.reindex(index=model.tag_names_, columns=model.time_index_) + z_loadings = loading[~is_traj] + z_loadings.index = z_loadings.index.get_level_values("tag") + return grid, z_loadings + + +def time_varying_loading_plot( + model: BatchPCA, + component: int = 1, + fig: go.Figure | None = None, + show_initial_conditions: bool = True, +) -> go.Figure: + """Plot one component's loadings as a function of time, one trace per tag. + + The batchwise-unfolded model has a separate loading for every + ``(tag, time)`` cell, so a component's loadings can be read as a set of + time-varying weight profiles: how strongly each variable loads on the + component at each point in the batch. Initial-condition (Z) loadings, which + have no time axis, are drawn as a marker group to the left of time zero. + + Parameters + ---------- + model : BatchPCA + A fitted :class:`process_improve.batch.BatchPCA` model. + component : int, default=1 + 1-based component index whose loadings to plot. + fig : plotly.graph_objects.Figure, optional + Figure to draw into; a new one is created when omitted. + show_initial_conditions : bool, default=True + Draw the initial-condition loadings (if the model has any) as a marker + group before time zero. + + Returns + ------- + plotly.graph_objects.Figure + """ + if not 0 < component <= model.n_components: + raise ValueError(f"The model has {model.n_components} components; need 1 <= component <= {model.n_components}.") + grid, z_loadings = _split_loadings(model, component) + + if fig is None: + fig = go.Figure() + for tag in model.tag_names_: + fig.add_trace( + go.Scatter( + x=list(model.time_index_), + y=grid.loc[tag].to_numpy(), + mode="lines", + name=str(tag), + ) + ) + if show_initial_conditions and len(z_loadings) > 0: + fig.add_trace( + go.Scatter( + x=[-1] * len(z_loadings), + y=z_loadings.to_numpy(), + mode="markers", + marker={"symbol": "diamond", "size": 9}, + text=[str(name) for name in z_loadings.index], + name="initial conditions", + ) + ) + fig.add_hline(y=0, line_color=REFERENCE_LINE_COLOR, line_width=1) + fig.update_layout( + template=DEFAULT_THEME, + title=f"Time-varying loadings for component {component}", + xaxis_title="Time [sequence order]", + yaxis_title=f"Loading p{component}", + ) + return fig + + +def contribution_at_time_plot( + contributions: pd.DataFrame, + k: int, + batch_id: typing.Hashable | None = None, + fig: go.Figure | None = None, +) -> go.Figure: + """Bar chart of per-tag contributions at one time sample, for one batch. + + Takes the output of :meth:`process_improve.batch.BatchPCA.spe_contributions` + or :meth:`~process_improve.batch.BatchPCA.t2_contributions` (one row per + batch, columns indexed by the 2-level ``(tag, sequence)`` unfolded index) + and shows, for a single batch and a single time sample ``k``, how much each + tag contributes. This localizes an abnormal event to the responsible + variable(s). + + Parameters + ---------- + contributions : pd.DataFrame + Contribution matrix from ``BatchPCA.spe_contributions`` / + ``t2_contributions``: one row per batch, a 2-level ``(tag, sequence)`` + column index. + k : int + The time sample (sequence value) at which to show the contributions. + batch_id : Hashable, optional + Which batch (row) to plot. Defaults to the first row; required to be a + valid row label when the matrix has more than one batch. + fig : plotly.graph_objects.Figure, optional + Figure to draw into; a new one is created when omitted. + + Returns + ------- + plotly.graph_objects.Figure + """ + if contributions.columns.nlevels != 2 or set(contributions.columns.names) != {"tag", "sequence"}: + raise ValueError( + "contributions must have a 2-level (tag, sequence) column index, as returned by " + "BatchPCA.spe_contributions / t2_contributions." + ) + if batch_id is None: + batch_id = contributions.index[0] + elif batch_id not in contributions.index: + raise ValueError(f"batch_id {batch_id!r} is not a row of the contributions matrix.") + + position = typing.cast("int", contributions.index.get_loc(batch_id)) + row = typing.cast("pd.Series", contributions.iloc[position]) + at_k = row[row.index.get_level_values("sequence") == k] + if at_k.empty: + raise ValueError(f"No contributions at time sample k={k}; available samples run over the batch length.") + tags = [str(label[0]) for label in at_k.index] + values = at_k.to_numpy(dtype=float) + + if fig is None: + fig = go.Figure() + fig.add_trace(go.Bar(x=tags, y=values, marker_color=np.where(values >= 0, "#2563EB", "#DC2626"))) + fig.add_hline(y=0, line_color=REFERENCE_LINE_COLOR, line_width=1) + fig.update_layout( + template=DEFAULT_THEME, + title=f"Contributions at time {k} (batch {batch_id})", + xaxis_title="Tag", + yaxis_title="Contribution", + ) + return fig diff --git a/tests/batch/test_batch_plots.py b/tests/batch/test_batch_plots.py new file mode 100644 index 00000000..1fe37f5e --- /dev/null +++ b/tests/batch/test_batch_plots.py @@ -0,0 +1,76 @@ +"""Tests for the batch model plots (time-varying loadings, contribution-at-time).""" + +import numpy as np +import pandas as pd +import pytest + +from process_improve.batch._batch_pca import BatchPCA +from process_improve.batch._batch_plots import contribution_at_time_plot, time_varying_loading_plot +from process_improve.batch.datasets import load_nylon +from process_improve.batch.preprocessing import resample_to_reference + +go = pytest.importorskip("plotly.graph_objects") + + +@pytest.fixture +def fitted_model() -> BatchPCA: + """Return a BatchPCA fitted on aligned nylon data.""" + batches = load_nylon() + tags = list(next(iter(batches.values())).columns) + aligned = resample_to_reference(batches, columns_to_align=tags, reference_batch=1) + return BatchPCA(n_components=3).fit(aligned) + + +def test_time_varying_loading_plot(fitted_model: BatchPCA) -> None: + """One trace per tag, each spanning the full batch length.""" + fig = time_varying_loading_plot(fitted_model, component=1) + assert isinstance(fig, go.Figure) + assert len(fig.data) == fitted_model.n_tags_ + assert len(fig.data[0].x) == fitted_model.n_timesteps_ + + +def test_time_varying_loading_plot_bad_component(fitted_model: BatchPCA) -> None: + """An out-of-range component index is rejected.""" + with pytest.raises(ValueError, match="components"): + time_varying_loading_plot(fitted_model, component=99) + + +def test_time_varying_loading_plot_with_initial_conditions() -> None: + """A model with a Z block adds a marker group for the initial conditions.""" + batches = load_nylon() + tags = list(next(iter(batches.values())).columns) + aligned = resample_to_reference(batches, columns_to_align=tags, reference_batch=1) + z = pd.DataFrame({"charge": [float(i) for i in range(len(aligned))]}, index=list(aligned.keys())) + model = BatchPCA(n_components=2).fit(aligned, initial_conditions=z) + fig = time_varying_loading_plot(model, component=1) + assert len(fig.data) == model.n_tags_ + 1 + + +def test_contribution_at_time_plot(fitted_model: BatchPCA) -> None: + """One bar per tag at the requested time sample.""" + batches = load_nylon() + tags = list(next(iter(batches.values())).columns) + aligned = resample_to_reference(batches, columns_to_align=tags, reference_batch=1) + scaled = fitted_model._scaled_wide(aligned, None) + contributions = fitted_model.spe_contributions(scaled) + fig = contribution_at_time_plot(contributions, k=57, batch_id=49) + assert isinstance(fig, go.Figure) + assert len(fig.data[0].x) == fitted_model.n_tags_ + + +def test_contribution_at_time_plot_defaults_to_first_batch(fitted_model: BatchPCA) -> None: + """With no batch_id the first row is used.""" + batches = load_nylon() + tags = list(next(iter(batches.values())).columns) + aligned = resample_to_reference(batches, columns_to_align=tags, reference_batch=1) + scaled = fitted_model._scaled_wide(aligned, None) + contributions = fitted_model.spe_contributions(scaled) + fig = contribution_at_time_plot(contributions, k=0) + assert len(fig.data[0].x) == fitted_model.n_tags_ + + +def test_contribution_at_time_plot_rejects_flat_columns() -> None: + """A contribution frame without the 2-level column index is rejected.""" + flat = pd.DataFrame(np.zeros((2, 3)), columns=["a", "b", "c"]) + with pytest.raises(ValueError, match="2-level"): + contribution_at_time_plot(flat, k=0) From 2b40e03d3d13bebe0d3c670363173994c297b172 Mon Sep 17 00:00:00 2001 From: Claude Date: Thu, 16 Jul 2026 21:23:43 +0000 Subject: [PATCH 7/8] Bump to 1.55.0: batch modelling, loaders, converters, changelog and docs Version bump for the new batch capabilities (BatchPCA, batch model plots, dataset loaders, completed data_input converters). CITATION.cff version and date-released synced to 1.55.0 / 2026-07-16; CHANGELOG section added with the compare-link footer; docs/api/batch.rst now autodocuments BatchPCA, the dataset loaders, and the batch plots. Co-Authored-By: Claude Fable 5 Claude-Session: https://claude.ai/code/session_01SyEpexCyqHn1rwTMYQSiP6 --- CHANGELOG.md | 34 +++++++++++++++++++++++++++++++++- CITATION.cff | 4 ++-- docs/api/batch.rst | 12 ++++++++++++ pyproject.toml | 2 +- 4 files changed, 48 insertions(+), 4 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 31ccf17d..10560f96 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -11,6 +11,37 @@ those changes. ## [Unreleased] +## [1.55.0] - 2026-07-16 + +### Added + +- `batch.BatchPCA`: batchwise-unfolded (multiway) PCA for batch trajectory + data. Unfolds an aligned batch-data dictionary batchwise (one row per batch), + optionally joins an initial-conditions (Z) block onto that row, mean-centres + and scales every column, and fits the existing `multivariate.PCA` by + composition. Exposes batch-level scores, SPE and Hotelling's T2 with control + limits, plus contribution and score/SPE/loading plot forwarders. This is the + Nomikos-MacGregor batch monitoring model. +- `batch.time_varying_loading_plot` and `batch.contribution_at_time_plot`: two + plotly plots for the unfolded model, showing a component's loadings over time + (one trace per tag) and the per-tag contribution to SPE or T2 at a chosen + time sample. +- `batch.load_nylon`, `batch.load_dryer`, and `batch.load_batch_fake_data`: + loader functions returning the bundled batch datasets as the standard + batch-data dictionary. + +### Changed + +- The `batch.data_input` converters `melted_to_wide`, `wide_to_dict`, and + `wide_to_melted` are now implemented (they previously returned an empty + dict / frame / `None`). They compose the existing conversions and round-trip + with `dict_to_wide` / `dict_to_melted`. + +### Fixed + +- `batch.plotting.plot_to_HTML` now passes the correct `responsive` Plotly + config key (previously the misspelled `resonsive` was silently ignored). + ## [1.54.0] - 2026-07-13 ### Changed @@ -2509,7 +2540,8 @@ this entry records them together. - Reworked the README with a sharper value proposition and a "Why not scikit-learn?" comparison table. -[Unreleased]: https://github.com/kgdunn/process-improve/compare/v1.54.0...HEAD +[Unreleased]: https://github.com/kgdunn/process-improve/compare/v1.55.0...HEAD +[1.55.0]: https://github.com/kgdunn/process-improve/compare/v1.54.0...v1.55.0 [1.54.0]: https://github.com/kgdunn/process-improve/compare/v1.53.0...v1.54.0 [1.53.0]: https://github.com/kgdunn/process-improve/compare/v1.52.4...v1.53.0 [1.52.4]: https://github.com/kgdunn/process-improve/compare/v1.52.3...v1.52.4 diff --git a/CITATION.cff b/CITATION.cff index 42e08d4b..e9a06570 100644 --- a/CITATION.cff +++ b/CITATION.cff @@ -12,8 +12,8 @@ authors: repository-code: "https://github.com/kgdunn/process-improve" url: "https://kgdunn.github.io/process-improve/" license: MIT -version: 1.54.0 -date-released: "2026-07-13" +version: 1.55.0 +date-released: "2026-07-16" keywords: - chemometrics - multivariate analysis diff --git a/docs/api/batch.rst b/docs/api/batch.rst index 1e7dd000..932b7969 100644 --- a/docs/api/batch.rst +++ b/docs/api/batch.rst @@ -1,6 +1,18 @@ Batch Data Analysis =================== +.. automodule:: process_improve.batch.datasets + :members: + :show-inheritance: + +.. autoclass:: process_improve.batch.BatchPCA + :members: + :show-inheritance: + +.. automodule:: process_improve.batch._batch_plots + :members: + :show-inheritance: + .. automodule:: process_improve.batch.features :members: :undoc-members: diff --git a/pyproject.toml b/pyproject.toml index cc772423..e7686464 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,6 +1,6 @@ [project] name = "process-improve" -version = "1.54.0" +version = "1.55.0" description = 'Designed Experiments; Latent Variables (PCA, PLS, multivariate methods with missing data); Process Monitoring; Batch data analysis.' readme = "README.md" license = "MIT" From 945710d7f640842a455aa4935f6f1e90802fa30e Mon Sep 17 00:00:00 2001 From: Claude Date: Thu, 16 Jul 2026 21:57:16 +0000 Subject: [PATCH 8/8] Cover BatchPCA and batch-plot error branches Adds tests for the validation and error paths flagged as uncovered: initial-conditions type/numeric/missing-value guards, group_by_batch with a Z block, fit_transform, diagnose with initial conditions, ellipse_coordinates, and the contribution-plot unknown-batch and out-of-range-time errors. Raises patch coverage of the new batch modules (_batch_pca to 97%, _batch_plots to 94%). Co-Authored-By: Claude Fable 5 Claude-Session: https://claude.ai/code/session_01SyEpexCyqHn1rwTMYQSiP6 --- tests/batch/test_batch_pca.py | 60 +++++++++++++++++++++++++++++++++ tests/batch/test_batch_plots.py | 20 +++++++++++ 2 files changed, 80 insertions(+) diff --git a/tests/batch/test_batch_pca.py b/tests/batch/test_batch_pca.py index 2f0290ff..cd564fa9 100644 --- a/tests/batch/test_batch_pca.py +++ b/tests/batch/test_batch_pca.py @@ -118,6 +118,66 @@ def test_transform_wrong_length_raises(aligned_nylon: dict) -> None: model.transform(shorter) +def test_initial_conditions_must_be_dataframe() -> None: + """A non-DataFrame Z block is rejected.""" + good = _synthetic_good_batches(n_batches=4) + with pytest.raises(TypeError, match="pandas DataFrame"): + BatchPCA(n_components=2).fit(good, initial_conditions={"charge": [1, 2, 3, 4]}) + + +def test_initial_conditions_must_be_numeric() -> None: + """A non-numeric Z column is rejected.""" + good = _synthetic_good_batches(n_batches=4) + z = pd.DataFrame({"grade": ["A", "B", "A", "C"]}, index=list(good.keys())) + with pytest.raises(ValueError, match="numeric"): + BatchPCA(n_components=2).fit(good, initial_conditions=z) + + +def test_initial_conditions_no_missing() -> None: + """A NaN in the Z block is rejected.""" + good = _synthetic_good_batches(n_batches=4) + z = pd.DataFrame({"charge": [1.0, np.nan, 3.0, 4.0]}, index=list(good.keys())) + with pytest.raises(ValueError, match="missing values"): + BatchPCA(n_components=2).fit(good, initial_conditions=z) + + +def test_group_by_batch_with_initial_conditions() -> None: + """group_by_batch labels Z loadings on the other level and still round-trips.""" + good = _synthetic_good_batches(n_batches=6) + z = pd.DataFrame({"charge": [float(i) for i in range(6)]}, index=list(good.keys())) + model = BatchPCA(n_components=2, group_by_batch=True).fit(good, initial_conditions=z) + assert model.n_initial_conditions_ == 1 + # with group_by_batch the Z label sits on the outer ("sequence") level as "" + z_rows = [label for label in model.loadings_.index if label[0] == ""] + assert {label[1] for label in z_rows} == {"charge"} + + +def test_fit_transform_returns_training_scores() -> None: + """fit_transform returns the same scores as fit followed by the attribute.""" + good = _synthetic_good_batches() + scores = BatchPCA(n_components=2).fit_transform(good) + reference = BatchPCA(n_components=2).fit(good).scores_ + np.testing.assert_allclose(scores.to_numpy(), reference.to_numpy(), atol=1e-8) + + +def test_diagnose_with_initial_conditions() -> None: + """Diagnosing accepts a Z block for the new batches.""" + good = _synthetic_good_batches() + z = pd.DataFrame({"charge": [float(i) for i in range(12)]}, index=list(good.keys())) + model = BatchPCA(n_components=2).fit(good, initial_conditions=z) + one = {"g0": good["g0"]} + z_new = pd.DataFrame({"charge": [0.0]}, index=["g0"]) + result = model.diagnose(one, initial_conditions=z_new) + assert result.scores.shape == (1, 2) + + +def test_ellipse_coordinates(aligned_nylon: dict) -> None: + """The T2 confidence ellipse returns paired coordinate arrays.""" + model = BatchPCA(n_components=3).fit(aligned_nylon) + x, y = model.ellipse_coordinates(1, 2) + assert len(x) == len(y) > 0 + + def test_unscaled_keeps_centring_only() -> None: """scale=False leaves the scaling factors at 1.0 (centring still applied).""" good = _synthetic_good_batches() diff --git a/tests/batch/test_batch_plots.py b/tests/batch/test_batch_plots.py index 1fe37f5e..cace198d 100644 --- a/tests/batch/test_batch_plots.py +++ b/tests/batch/test_batch_plots.py @@ -74,3 +74,23 @@ def test_contribution_at_time_plot_rejects_flat_columns() -> None: flat = pd.DataFrame(np.zeros((2, 3)), columns=["a", "b", "c"]) with pytest.raises(ValueError, match="2-level"): contribution_at_time_plot(flat, k=0) + + +def test_contribution_at_time_plot_unknown_batch_id(fitted_model: BatchPCA) -> None: + """An unknown batch_id is rejected with a clear error.""" + batches = load_nylon() + tags = list(next(iter(batches.values())).columns) + aligned = resample_to_reference(batches, columns_to_align=tags, reference_batch=1) + contributions = fitted_model.spe_contributions(fitted_model._scaled_wide(aligned, None)) + with pytest.raises(ValueError, match="not a row"): + contribution_at_time_plot(contributions, k=0, batch_id=99999) + + +def test_contribution_at_time_plot_out_of_range_time(fitted_model: BatchPCA) -> None: + """A time sample beyond the batch length has no contributions and is rejected.""" + batches = load_nylon() + tags = list(next(iter(batches.values())).columns) + aligned = resample_to_reference(batches, columns_to_align=tags, reference_batch=1) + contributions = fitted_model.spe_contributions(fitted_model._scaled_wide(aligned, None)) + with pytest.raises(ValueError, match="No contributions at time"): + contribution_at_time_plot(contributions, k=fitted_model.n_timesteps_ + 100)