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" diff --git a/src/process_improve/batch/__init__.py b/src/process_improve/batch/__init__.py index 5bb3e272..4601862a 100644 --- a/src/process_improve/batch/__init__.py +++ b/src/process_improve/batch/__init__.py @@ -1,13 +1,24 @@ -"""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._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, 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, @@ -36,10 +47,13 @@ ) __all__ = [ + # Modelling + "BatchPCA", # Preprocessing/alignment "batch_dtw", # Data input/conversion "check_valid_batch_dict", + "contribution_at_time_plot", "cross", "determine_scaling", "dict_to_melted", @@ -63,8 +77,14 @@ "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", + "time_varying_loading_plot", + "wide_to_dict", "wide_to_melted", ] 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/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/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/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/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( diff --git a/tests/batch/test_batch_pca.py b/tests/batch/test_batch_pca.py new file mode 100644 index 00000000..cd564fa9 --- /dev/null +++ b/tests/batch/test_batch_pca.py @@ -0,0 +1,187 @@ +"""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_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() + 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) diff --git a/tests/batch/test_batch_plots.py b/tests/batch/test_batch_plots.py new file mode 100644 index 00000000..cace198d --- /dev/null +++ b/tests/batch/test_batch_plots.py @@ -0,0 +1,96 @@ +"""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) + + +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) 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: 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