Utilities

Standardization, statistics helpers, and the arviz/pymc version-drift compatibility shims (route raw arviz calls through arviz_compat).

Utility modules for MMM Framework.

class mmm_framework.utils.DataStandardizer(epsilon=1e-08)[source]

Bases: object

Standardize data with zero mean and unit variance.

This class provides methods for standardizing data (z-score normalization) which is essential for Bayesian models. It handles both 1D and 2D data, and includes a small epsilon term to prevent division by zero for constant data.

Parameters

epsilonfloat, default=1e-8

Small constant added to standard deviation to prevent division by zero.

Examples

>>> import numpy as np
>>> from mmm_framework.utils import DataStandardizer
>>>
>>> # Create standardizer
>>> standardizer = DataStandardizer()
>>>
>>> # Fit and transform training data
>>> data = np.array([10.0, 20.0, 30.0, 40.0, 50.0])
>>> standardized, params = standardizer.fit_transform(data)
>>>
>>> # Transform new data using same parameters
>>> new_data = np.array([25.0, 35.0])
>>> transformed = standardizer.transform(new_data, params)
>>>
>>> # Reverse transformation
>>> original_scale = standardizer.inverse_transform(transformed, params)
__init__(epsilon=1e-08)[source]

Initialize DataStandardizer.

Parameters

epsilonfloat, default=1e-8

Small constant added to standard deviation to prevent division by zero when data has zero variance.

fit(data)[source]

Compute standardization parameters from data.

Return type:

StandardizationParams

Parameters

dataNDArray

Input data to compute parameters from. Can be 1D or 2D. For 2D data, parameters are computed per column (axis=0).

Returns

StandardizationParams

Parameters containing mean and standard deviation.

transform(data, params=None)[source]

Apply standardization to data.

Return type:

numpy.ndarray

Parameters

dataNDArray

Data to standardize.

paramsStandardizationParams, optional

Parameters to use for transformation. If None, uses parameters from most recent fit() call.

Returns

NDArray

Standardized data with zero mean and unit variance.

Raises

ValueError

If params is None and fit() has not been called.

fit_transform(data)[source]

Fit parameters and transform data in one step.

Return type:

tuple[numpy.ndarray, StandardizationParams]

Parameters

dataNDArray

Data to fit and transform.

Returns

tuple[NDArray, StandardizationParams]

Tuple of (standardized_data, parameters).

inverse_transform(data, params=None)[source]

Reverse standardization to recover original scale.

Return type:

numpy.ndarray

Parameters

dataNDArray

Standardized data to transform back.

paramsStandardizationParams, optional

Parameters to use for inverse transformation. If None, uses parameters from most recent fit() call.

Returns

NDArray

Data in original scale.

Raises

ValueError

If params is None and fit() has not been called.

class mmm_framework.utils.StandardizationParams(mean, std)[source]

Bases: object

Parameters from standardization fit.

Stores the mean and standard deviation used for standardization, allowing the transformation to be applied to new data or reversed.

Attributes

meanfloat | NDArray

Mean value(s) used for centering. Scalar for 1D data, array for multi-dimensional data.

stdfloat | NDArray

Standard deviation(s) used for scaling. Scalar for 1D data, array for multi-dimensional data.

mean: float | NDArray
std: float | NDArray
to_dict()[source]

Convert to serializable dictionary.

Return type:

dict

Returns

dict

Dictionary with ‘mean’ and ‘std’ keys, with numpy arrays converted to lists for JSON serialization.

classmethod from_dict(d)[source]

Create from dictionary.

Return type:

StandardizationParams

Parameters

ddict

Dictionary with ‘mean’ and ‘std’ keys.

Returns

StandardizationParams

Reconstructed parameters object.

__init__(mean, std)
mmm_framework.utils.compute_hdi_bounds(samples, hdi_prob=0.94, axis=0)[source]

Compute highest density interval bounds using percentiles.

Computes the central credible interval bounds for a given probability mass. This uses a simple percentile-based approach which is appropriate for approximately symmetric distributions.

Return type:

tuple[numpy.ndarray, numpy.ndarray]

Parameters

samplesNDArray

Sample array from posterior distribution. Shape can be (n_samples,) for 1D or (n_samples, n_observations) for 2D.

hdi_probfloat, default=0.94

Probability mass for the HDI. For example, 0.94 gives the central 94% interval.

axisint, default=0

Axis along which to compute percentiles. Typically axis=0 when samples are in the first dimension.

Returns

tuple[NDArray, NDArray]

Tuple of (lower_bound, upper_bound) arrays. Shape depends on input shape and axis parameter.

Examples

>>> import numpy as np
>>> from mmm_framework.utils import compute_hdi_bounds
>>>
>>> # Generate samples
>>> np.random.seed(42)
>>> samples = np.random.randn(1000, 10)  # 1000 samples, 10 observations
>>>
>>> # Compute 94% HDI
>>> lower, upper = compute_hdi_bounds(samples, hdi_prob=0.94)
>>> print(f"Lower bounds shape: {lower.shape}")  # (10,)
>>> print(f"Upper bounds shape: {upper.shape}")  # (10,)

Notes

This function uses a simple percentile-based approach rather than a true highest density interval algorithm. For symmetric distributions like the Normal distribution, this is equivalent to the HDI. For highly skewed distributions, a proper HDI algorithm may give different (narrower) intervals.

The percentiles are computed as: - lower = (1 - hdi_prob) / 2 * 100 - upper = (1 + hdi_prob) / 2 * 100

For hdi_prob=0.94, this gives percentiles 3 and 97.

mmm_framework.utils.logged_suppress(context, *exceptions, level='DEBUG')[source]

Suppress exceptions (default: Exception) but log them with context.

Parameters

contextstr

Short description of the optional work being attempted, included in the log line so a swallowed failure can be traced.

*exceptionsexception types

Which exceptions to suppress. Defaults to Exception (not BaseException — KeyboardInterrupt/SystemExit always propagate).

levelstr

loguru level for the log line (default DEBUG).

Arviz Compat

Version-robust shims for arviz / pymc API drift.

The framework targets a range of arviz and pymc versions whose public APIs shifted in ways that fail loudly (TypeError/AttributeError) or, worse, silently (wrong result). Centralizing the shims here keeps every call site honest and avoids re-discovering the same breakage module by module.

Covered drift: - pm.sample_prior_predictive renamed samples -> draws (pymc >=5.x). - arviz containers migrated InferenceData -> xarray DataTree (arviz

>=0.22): .groups became a property of slash-prefixed paths (not a method), .extend is gone, and Dataset.to_array became to_dataarray.

  • az.from_dict flipped its calling convention across that migration AND the wrong form fails silently (wraps everything as one var "posterior").

  • az.hdi renamed hdi_prob -> prob (arviz >=1.x) AND changed 2-d ndarray semantics from pooled (chain, draw) to a batch of independent rows — use hdi_bounds() for the historical pooled interval.

mmm_framework.utils.arviz_compat.sample_prior_predictive(samples, random_seed=None)[source]

Call pm.sample_prior_predictive across the samples->``draws`` rename.

Must be called inside a with model: context. Callers keep their public samples argument.

mmm_framework.utils.arviz_compat.dataset_extremum(ds, kind)[source]

Reduce an arviz per-variable stats container to one scalar.

az.rhat / az.ess return an xarray Dataset on legacy arviz and a DataTree on arviz >=0.22; both expose data_vars, so we reduce over those rather than relying on the renamed to_array/to_dataarray. All-NaN variables (e.g. R-hat for a single-chain fit) are skipped so the reduction stays quiet and meaningful. kind is "max" or "min".

Return type:

float

mmm_framework.utils.arviz_compat.hdi_bounds(samples, prob)[source]

(lo, hi) highest-density interval of the POOLED samples.

Absorbs two arviz 1.x changes at once: hdi_prob was renamed prob, and a 2-d ndarray input is now treated as a BATCH of independent rows (one interval per row) instead of pooled (chain, draw) draws. Inputs are flattened here so every call site keeps the historical pooled semantics.

Return type:

tuple[float, float]

mmm_framework.utils.arviz_compat.group_names(idata)[source]

Return an arviz container’s group names, normalized (no leading /).

Legacy InferenceData exposes groups() as a method returning bare names; the newer DataTree exposes groups as a property of slash-prefixed paths. Both are normalized to bare names here.

Return type:

list[str]

mmm_framework.utils.arviz_compat.has_group(idata, name)[source]

True if an arviz container exposes the named group (robust across APIs).

Return type:

bool

mmm_framework.utils.arviz_compat.attach_prior(trace, prior)[source]

Merge the prior groups into the posterior trace, best-effort.

InferenceData.extend exists on the legacy container but not on the newer DataTree. Prior groups are only used for prior-vs-posterior tooling, so on any incompatibility we warn and return the trace unchanged.

mmm_framework.utils.arviz_compat.posterior_from_dict(posterior)[source]

Build an arviz container whose posterior group holds posterior.

Values must already be shaped (chain, draw, *shape). arviz flipped from_dict’s calling convention between the legacy InferenceData and the DataTree era — and the wrong form does not raise: it silently wraps everything as one variable named "posterior". So we try both conventions and keep whichever actually materialized the variables; if neither does, build the dataset by hand.

mmm_framework.utils.arviz_compat.summary(data, var_names=None, **kwargs)[source]

az.summary with NUMERIC values across the arviz 1.x formatting change.

arviz 1.x defaults round_to="auto", which returns formatted STRINGS (object dtype) — silent poison for numeric consumers: max() on the r_hat column becomes a LEXICOGRAPHIC string max (“9.99” > “10.01”). round_to="none" restores raw floats on both arviz lines. Note the interval columns also drifted (hdi_3%/hdi_97%eti89_lb/ eti89_ub, an 89% equal-tailed default) — pass ci_prob/ci_kind explicitly if you consume them.

mmm_framework.utils.arviz_compat.hdi_dataset(data, prob, var_names=None)[source]

Trace-level HDI across the hdi_probprob rename.

result[var].values keeps the legacy contract on both arviz lines: shape (*var_shape, 2) with the last axis = (lower, upper).

mmm_framework.utils.arviz_compat.plot_posterior(data, var_names=None, **kwargs)[source]

Posterior-distribution plot across the arviz-plots split.

Legacy arviz exposed az.plot_posterior; arviz 1.x moved plotting into arviz_plots where the equivalent is plot_dist (defaults to the posterior group). az.plot_trace survived the split (re-exported), so only this one needs a shim.

mmm_framework.utils.arviz_compat.dataset_to_idata(ds, group='posterior')[source]

Wrap a ready xarray Dataset as one group of an arviz container.

az.InferenceData(posterior=ds) stopped working on the DataTree migration (DataTree.__init__ takes no group kwargs); the DataTree era assigns groups by key instead. Works on both.

mmm_framework.utils.arviz_compat.point_to_idata(point)[source]

Wrap a find_MAP point dict into a (chain=1, draw=1) InferenceData.

Drops PyTensor’s transformed duplicates (*_log__ / *_interval__) and keeps the constrained values plus deterministics, matching the variable names produced by NUTS. Container construction (and the from_dict convention flip) is handled by posterior_from_dict().

mmm_framework.utils.arviz_compat.psis_log_weights(log_lik)[source]

PSIS-smoothed LOO log-weights from pointwise log-likelihood draws.

log_lik is (n_obs, n_samples) (draws on the last axis). Returns (log_weights, khat) with log_weights the same shape (normalized so logsumexp == 0 per observation) and khat of shape (n_obs,) — the Pareto shape diagnostic (khat > 0.7 means the weights for that observation are unreliable).

arviz 1.x removed the top-level az.psislw; the smoother now lives on arviz_stats.base.array_stats. The legacy az.psislw (draws on the FIRST axis) is the fallback for pre-1.x environments.

Return type:

tuple[ndarray, ndarray]

Standardization

Data standardization utilities for MMM Framework.

This module provides utilities for standardizing data (zero mean, unit variance) which is commonly needed for Bayesian models to ensure numerical stability.

class mmm_framework.utils.standardization.StandardizationParams(mean, std)[source]

Bases: object

Parameters from standardization fit.

Stores the mean and standard deviation used for standardization, allowing the transformation to be applied to new data or reversed.

Attributes

meanfloat | NDArray

Mean value(s) used for centering. Scalar for 1D data, array for multi-dimensional data.

stdfloat | NDArray

Standard deviation(s) used for scaling. Scalar for 1D data, array for multi-dimensional data.

mean: float | NDArray
std: float | NDArray
to_dict()[source]

Convert to serializable dictionary.

Return type:

dict

Returns

dict

Dictionary with ‘mean’ and ‘std’ keys, with numpy arrays converted to lists for JSON serialization.

classmethod from_dict(d)[source]

Create from dictionary.

Return type:

StandardizationParams

Parameters

ddict

Dictionary with ‘mean’ and ‘std’ keys.

Returns

StandardizationParams

Reconstructed parameters object.

__init__(mean, std)
class mmm_framework.utils.standardization.DataStandardizer(epsilon=1e-08)[source]

Bases: object

Standardize data with zero mean and unit variance.

This class provides methods for standardizing data (z-score normalization) which is essential for Bayesian models. It handles both 1D and 2D data, and includes a small epsilon term to prevent division by zero for constant data.

Parameters

epsilonfloat, default=1e-8

Small constant added to standard deviation to prevent division by zero.

Examples

>>> import numpy as np
>>> from mmm_framework.utils import DataStandardizer
>>>
>>> # Create standardizer
>>> standardizer = DataStandardizer()
>>>
>>> # Fit and transform training data
>>> data = np.array([10.0, 20.0, 30.0, 40.0, 50.0])
>>> standardized, params = standardizer.fit_transform(data)
>>>
>>> # Transform new data using same parameters
>>> new_data = np.array([25.0, 35.0])
>>> transformed = standardizer.transform(new_data, params)
>>>
>>> # Reverse transformation
>>> original_scale = standardizer.inverse_transform(transformed, params)
__init__(epsilon=1e-08)[source]

Initialize DataStandardizer.

Parameters

epsilonfloat, default=1e-8

Small constant added to standard deviation to prevent division by zero when data has zero variance.

fit(data)[source]

Compute standardization parameters from data.

Return type:

StandardizationParams

Parameters

dataNDArray

Input data to compute parameters from. Can be 1D or 2D. For 2D data, parameters are computed per column (axis=0).

Returns

StandardizationParams

Parameters containing mean and standard deviation.

transform(data, params=None)[source]

Apply standardization to data.

Return type:

numpy.ndarray

Parameters

dataNDArray

Data to standardize.

paramsStandardizationParams, optional

Parameters to use for transformation. If None, uses parameters from most recent fit() call.

Returns

NDArray

Standardized data with zero mean and unit variance.

Raises

ValueError

If params is None and fit() has not been called.

fit_transform(data)[source]

Fit parameters and transform data in one step.

Return type:

tuple[numpy.ndarray, StandardizationParams]

Parameters

dataNDArray

Data to fit and transform.

Returns

tuple[NDArray, StandardizationParams]

Tuple of (standardized_data, parameters).

inverse_transform(data, params=None)[source]

Reverse standardization to recover original scale.

Return type:

numpy.ndarray

Parameters

dataNDArray

Standardized data to transform back.

paramsStandardizationParams, optional

Parameters to use for inverse transformation. If None, uses parameters from most recent fit() call.

Returns

NDArray

Data in original scale.

Raises

ValueError

If params is None and fit() has not been called.

Statistics

Statistical utility functions for MMM Framework.

This module provides statistical utilities commonly used in Bayesian model analysis, such as computing highest density intervals (HDI).

mmm_framework.utils.statistics.compute_hdi_bounds(samples, hdi_prob=0.94, axis=0)[source]

Compute highest density interval bounds using percentiles.

Computes the central credible interval bounds for a given probability mass. This uses a simple percentile-based approach which is appropriate for approximately symmetric distributions.

Return type:

tuple[numpy.ndarray, numpy.ndarray]

Parameters

samplesNDArray

Sample array from posterior distribution. Shape can be (n_samples,) for 1D or (n_samples, n_observations) for 2D.

hdi_probfloat, default=0.94

Probability mass for the HDI. For example, 0.94 gives the central 94% interval.

axisint, default=0

Axis along which to compute percentiles. Typically axis=0 when samples are in the first dimension.

Returns

tuple[NDArray, NDArray]

Tuple of (lower_bound, upper_bound) arrays. Shape depends on input shape and axis parameter.

Examples

>>> import numpy as np
>>> from mmm_framework.utils import compute_hdi_bounds
>>>
>>> # Generate samples
>>> np.random.seed(42)
>>> samples = np.random.randn(1000, 10)  # 1000 samples, 10 observations
>>>
>>> # Compute 94% HDI
>>> lower, upper = compute_hdi_bounds(samples, hdi_prob=0.94)
>>> print(f"Lower bounds shape: {lower.shape}")  # (10,)
>>> print(f"Upper bounds shape: {upper.shape}")  # (10,)

Notes

This function uses a simple percentile-based approach rather than a true highest density interval algorithm. For symmetric distributions like the Normal distribution, this is equivalent to the HDI. For highly skewed distributions, a proper HDI algorithm may give different (narrower) intervals.

The percentiles are computed as: - lower = (1 - hdi_prob) / 2 * 100 - upper = (1 + hdi_prob) / 2 * 100

For hdi_prob=0.94, this gives percentiles 3 and 97.