Design

Experimental design: factors, grids, and screening.

Sobol evaluation budgets

For base sample size N and D continuous factors, both screen(method="sobol", n_trajectories=N) and sobol_indices(n_samples=N) currently request the second-order Saltelli design with N * (2*D + 2) nominal configurations. screen returns first-order indices, but still uses that larger design. A completed run invokes run_fn N * (2*D + 2) * n_reps times. Repeated configurations in the design do not reduce that callback count.

For example, N=8, D=3 and n_reps=2 mean 64 nominal configurations and 128 callable invocations. Replication costs apply even when a callback does not accept the optional rep keyword; meaningful stochastic replication requires the caller's randomness contract. These counts describe the screening callback, not hidden work or caching inside it. Surrogate fitting/CV has its own cost. This correction preserves the sampler/analyzer settings and existing defaults; there is no new first/total-only sampling option.

trade_study.FactorType

Bases: Enum

Type of design factor.

trade_study.Factor(name, factor_type, levels=None, bounds=None, log_scale=False) dataclass

A single design factor.

Attributes:

Name	Type	Description
name	str	Factor identifier (e.g. "alpha", "layer1_method").
factor_type	FactorType	Continuous, discrete, or categorical.
levels	list[Any] | None	For categorical/discrete: list of allowed values.
bounds	tuple[float, float] | None	For continuous: (low, high) tuple.
log_scale	bool	For continuous factors with positive bounds: sample and model the factor uniformly in log space, for parameters that span orders of magnitude.

__post_init__()

Validate factor constraints.

Raises:

Type	Description
ValueError	If name is empty, continuous factor has missing or invalid bounds, a log-scale factor has non-positive bounds or is not continuous, or discrete/categorical factor has empty levels.

Source code in src/trade_study/design.py

def __post_init__(self) -> None:
    """Validate factor constraints.

    Raises:
        ValueError: If name is empty, continuous factor has missing or
            invalid bounds, a log-scale factor has non-positive bounds or
            is not continuous, or discrete/categorical factor has empty
            levels.
    """
    if not self.name:
        msg = "Factor name must be a non-empty string"
        raise ValueError(msg)
    if self.factor_type == FactorType.CONTINUOUS:
        if self.bounds is None:
            msg = f"Continuous factor '{self.name}' requires bounds"
            raise ValueError(msg)
        lo, hi = self.bounds
        if not (np.isfinite(lo) and np.isfinite(hi)):
            msg = f"Continuous factor '{self.name}' bounds must be finite"
            raise ValueError(msg)
        if lo >= hi:
            msg = f"Continuous factor '{self.name}' requires lo < hi"
            raise ValueError(msg)
        if self.log_scale and lo <= 0:
            msg = f"Log-scale factor '{self.name}' requires positive bounds"
            raise ValueError(msg)
    else:
        if self.log_scale:
            msg = f"Factor '{self.name}': log_scale applies to continuous factors"
            raise ValueError(msg)
        if self.levels is None:
            msg = f"Factor '{self.name}' of type {self.factor_type} requires levels"
            raise ValueError(msg)
        if len(self.levels) == 0:
            msg = f"Factor '{self.name}' levels must be non-empty"
            raise ValueError(msg)

trade_study.FactorConstraint(predicate, name='factor_constraint') dataclass

Coupled-factor constraint applied at design-generation time.

Unlike :class:trade_study.protocols.Constraint, which filters trials after simulation based on observables, a FactorConstraint filters candidate configurations before they are evaluated. It is the way to express coupled constraints such as "if method == 'a' then patience <= 2" without distorting the design via post-hoc filtering of returned results.

Attributes:

Name	Type	Description
predicate	Callable[[dict[str, Any]], bool]	Callable taking a candidate config dict and returning True when the config is feasible.
name	str	Human-readable label used in error messages and warnings.

is_feasible(config)

Return whether config satisfies this constraint.

Parameters:

Name	Type	Description	Default
config	dict[str, Any]	Candidate design point.	required

Returns:

Type	Description
bool	True if the predicate accepts the config.

Source code in src/trade_study/design.py

def is_feasible(self, config: dict[str, Any]) -> bool:
    """Return whether ``config`` satisfies this constraint.

    Args:
        config: Candidate design point.

    Returns:
        ``True`` if the predicate accepts the config.
    """
    return bool(self.predicate(config))

trade_study.build_grid(factors, *, method='full', n_samples=100, seed=42, scramble=True, constraints=None)

Build an experimental design grid.

Parameters:

Name	Type	Description	Default
factors	list[Factor]	List of design factors.	required
method	str	Design method. One of: - "full": Full factorial (categorical/discrete only). - "lhs": Latin hypercube sampling (continuous factors, maps categorical factors to uniform random selection). - "sobol": Scrambled Sobol' sequence via scipy.stats.qmc. - "halton": Scrambled Halton sequence via scipy.stats.qmc.	'full'
n_samples	int	Number of samples for LHS / QMC methods.	100
seed	int	Random seed.	42
scramble	bool	Whether to apply scrambling to QMC sequences (Sobol / Halton). Ignored for other methods.	True
constraints	Sequence[FactorConstraint] | None	Optional sequence of :class:FactorConstraint applied at design-generation time. For "full", the Cartesian product is filtered. For "lhs"/"sobol"/ "halton", rejection sampling with oversampling is used to return n_samples feasible configs (rejection breaks the QMC space-filling guarantee; expect higher discrepancy when the feasibility ratio is low).	None

Returns:

Type	Description
list[dict[str, Any]]	List of config dictionaries, one per design point.

Raises:

Type	Description
ValueError	If an unknown design method is specified, or if constraints reject more than the maximum oversample budget allows.

Examples:

>>> factors = [Factor("size", FactorType.DISCRETE, levels=[10, 20])]
>>> build_grid(factors)
[{'size': 10}, {'size': 20}]
>>> constraint = FactorConstraint(lambda c: c["size"] <= 10)
>>> build_grid(factors, constraints=[constraint])
[{'size': 10}]

Source code in src/trade_study/design.py

def build_grid(
    factors: list[Factor],
    *,
    method: str = "full",
    n_samples: int = 100,
    seed: int = 42,
    scramble: bool = True,
    constraints: Sequence[FactorConstraint] | None = None,
) -> list[dict[str, Any]]:
    """Build an experimental design grid.

    Args:
        factors: List of design factors.
        method: Design method. One of:
            - "full": Full factorial (categorical/discrete only).
            - "lhs": Latin hypercube sampling (continuous factors, maps
              categorical factors to uniform random selection).
            - "sobol": Scrambled Sobol' sequence via ``scipy.stats.qmc``.
            - "halton": Scrambled Halton sequence via ``scipy.stats.qmc``.
        n_samples: Number of samples for LHS / QMC methods.
        seed: Random seed.
        scramble: Whether to apply scrambling to QMC sequences (Sobol /
            Halton). Ignored for other methods.
        constraints: Optional sequence of :class:`FactorConstraint`
            applied at design-generation time. For ``"full"``, the
            Cartesian product is filtered. For ``"lhs"``/``"sobol"``/
            ``"halton"``, rejection sampling with oversampling is used
            to return ``n_samples`` feasible configs (rejection breaks
            the QMC space-filling guarantee; expect higher discrepancy
            when the feasibility ratio is low).

    Returns:
        List of config dictionaries, one per design point.

    Raises:
        ValueError: If an unknown design method is specified, or if
            constraints reject more than the maximum oversample budget
            allows.

    Examples:
        >>> factors = [Factor("size", FactorType.DISCRETE, levels=[10, 20])]
        >>> build_grid(factors)
        [{'size': 10}, {'size': 20}]
        >>> constraint = FactorConstraint(lambda c: c["size"] <= 10)
        >>> build_grid(factors, constraints=[constraint])
        [{'size': 10}]
    """
    if method == "full":
        return _full_factorial(factors, constraints=constraints)
    if method == "lhs":
        return _latin_hypercube(
            factors,
            n_samples=n_samples,
            seed=seed,
            constraints=constraints,
        )
    if method in {"sobol", "halton"}:
        return _qmc_sample(
            factors,
            n_samples=n_samples,
            seed=seed,
            qmc_method=method,
            scramble=scramble,
            constraints=constraints,
        )
    msg = f"Unknown design method: {method!r}"
    raise ValueError(msg)

trade_study.reduce_factors(factors, importance, *, threshold=0.1)

Keep only factors whose max importance exceeds threshold.

Parameters:

Name	Type	Description	Default
factors	list[Factor]	Original factor list.	required
importance	dict[str, NDArray[floating[Any]]]	Output of screen().	required
threshold	float	Minimum importance to retain a factor.	0.1

Returns:

Type	Description
list[Factor]	Reduced list of influential factors.

Warns:

Type	Description
UserWarning	If every observable's importance is NaN for one or more factors -- those factors are dropped (there's no valid data to compare against threshold), but that reflects missing data, not confirmed unimportance.

Source code in src/trade_study/design.py

def reduce_factors(
    factors: list[Factor],
    importance: dict[str, NDArray[np.floating[Any]]],
    *,
    threshold: float = 0.1,
) -> list[Factor]:
    """Keep only factors whose max importance exceeds threshold.

    Args:
        factors: Original factor list.
        importance: Output of ``screen()``.
        threshold: Minimum importance to retain a factor.

    Returns:
        Reduced list of influential factors.

    Warns:
        UserWarning: If every observable's importance is NaN for one or
            more factors -- those factors are dropped (there's no
            valid data to compare against ``threshold``), but that
            reflects missing data, not confirmed unimportance.
    """
    continuous = [f for f in factors if f.factor_type == FactorType.CONTINUOUS]
    non_continuous = [f for f in factors if f.factor_type != FactorType.CONTINUOUS]

    # NaN-filled, not zero-filled: a factor with zero valid (non-NaN)
    # measurements across every observable must stay NaN throughout, not
    # silently settle at 0.0 -- indistinguishable from "confirmed
    # unimportant" otherwise. np.fmax (not np.maximum) ignores NaN in
    # either operand rather than propagating it, so one observable that's
    # legitimately undefined in some regimes (e.g. a Type-I rate, NaN
    # outside null regimes) can't erase a real, significant importance
    # value a *different* observable found for the same factor (#119).
    max_importance = np.full(len(continuous), np.nan)
    for arr in importance.values():
        max_importance = np.fmax(max_importance, arr)

    all_nan = np.isnan(max_importance)
    if all_nan.any():
        names = [continuous[i].name for i in np.flatnonzero(all_nan)]
        warnings.warn(
            f"reduce_factors: every observable's importance is NaN for "
            f"{names} -- dropping them from the reduced list, but this "
            f"reflects missing data, not confirmed unimportance. See "
            f"https://github.com/jcm-sci/trade-study/issues/119.",
            UserWarning,
            stacklevel=2,
        )

    kept = [
        f for f, imp in zip(continuous, max_importance, strict=True) if imp >= threshold
    ]
    return non_continuous + kept

trade_study.screen(run_fn, factors, *, method='morris', n_trajectories=100, n_reps=1, seed=42)

Screen factors for influence on observables via SALib.

Parameters:

Name	Type	Description	Default
run_fn	Callable[[dict[str, Any]], dict[str, float]]	Callable that takes a config dict and returns a dict of observable name → scalar score.	required
factors	list[Factor]	List of continuous factors to screen.	required
method	str	Screening method ("morris" or "sobol").	'morris'
n_trajectories	int	Number of Morris trajectories. For Sobol, this controls the base sample size N; the total number of model configurations is N x (2 x num_vars + 2), using the current second-order Saltelli design even though screen returns S1 only. Physical run_fn calls multiply this by n_reps.	100
n_reps	int	Replicate draws averaged into each sampled point's scores before SALib sees them, the same convention run_grid(..., n_reps>1) uses -- run_fn opts in via an optional keyword-only rep parameter, detected via introspection. Without it, a single stochastic draw per point can produce importance estimates that are artifacts of that one draw rather than the factor's real effect. Default 1 (a single draw per point, today's behavior).	1
seed	int	Random seed.	42

Returns:

Type	Description
dict[str, NDArray[floating[Any]]]	Dictionary mapping observable names to arrays of factor importance
dict[str, NDArray[floating[Any]]]	(mu_star for Morris, S1 for Sobol), one value per factor.

Raises:

Type	Description
ValueError	If method is unknown, no continuous factors are provided, or n_reps is less than 1.

Source code in src/trade_study/design.py

def screen(
    run_fn: Callable[[dict[str, Any]], dict[str, float]],
    factors: list[Factor],
    *,
    method: str = "morris",
    n_trajectories: int = 100,
    n_reps: int = 1,
    seed: int = 42,
) -> dict[str, NDArray[np.floating[Any]]]:
    """Screen factors for influence on observables via SALib.

    Args:
        run_fn: Callable that takes a config dict and returns a dict of
            observable name → scalar score.
        factors: List of continuous factors to screen.
        method: Screening method (``"morris"`` or ``"sobol"``).
        n_trajectories: Number of Morris trajectories.  For Sobol, this
            controls the base sample size *N*; the total number of model
            configurations is *N* x (2 x num_vars + 2), using the current
            second-order Saltelli design even though ``screen`` returns S1
            only. Physical ``run_fn`` calls multiply this by ``n_reps``.
        n_reps: Replicate draws averaged into each sampled point's scores
            before SALib sees them, the same convention
            ``run_grid(..., n_reps>1)`` uses -- ``run_fn`` opts in
            via an optional keyword-only ``rep`` parameter, detected via
            introspection. Without it, a single stochastic draw per point
            can produce importance estimates that are artifacts of that
            one draw rather than the factor's real effect. Default 1 (a
            single draw per point, today's behavior).
        seed: Random seed.

    Returns:
        Dictionary mapping observable names to arrays of factor importance
        (mu_star for Morris, S1 for Sobol), one value per factor.

    Raises:
        ValueError: If *method* is unknown, no continuous factors are
            provided, or ``n_reps`` is less than 1.
    """
    continuous = [f for f in factors if f.factor_type == FactorType.CONTINUOUS]
    if not continuous:
        msg = "Screening requires at least one continuous factor"
        raise ValueError(msg)
    if n_reps < 1:
        msg = f"n_reps must be >= 1; got {n_reps}"
        raise ValueError(msg)

    problem = _salib_problem(continuous)

    if method == "morris":
        return _screen_morris(run_fn, problem, n_trajectories, n_reps, seed)
    if method == "sobol":
        return _screen_sobol(run_fn, problem, n_trajectories, n_reps, seed)

    msg = f"Unknown screening method: {method!r}"
    raise ValueError(msg)

trade_study.sobol_indices(run_fn, factors, *, n_samples=100, n_reps=1, seed=42)

Sobol first- and total-order sensitivity indices.

Like screen(method="sobol"), but returns both S1 (first-order) and ST (total-order) per observable instead of discarding ST. ST - S1 is the standard way to detect interaction effects: a factor with small S1 but large ST is interacting with other factors rather than acting independently -- information screen() can't surface at all, since it only keeps S1.

Parameters:

Name	Type	Description	Default
run_fn	Callable[[dict[str, Any]], dict[str, float]]	Callable that takes a config dict and returns a dict of observable name -> scalar score.	required
factors	list[Factor]	Factor list to analyze; only continuous factors are varied (as in screen()).	required
n_samples	int	Base sample size N for the second-order Saltelli design; total configurations are N x (2 x num_vars + 2). Physical run_fn calls multiply this by n_reps.	100
n_reps	int	Replicate draws averaged into each sampled point's scores before SALib sees them -- see screen()'s n_reps for the full rationale and the run_fn opt-in convention. Default 1 (a single draw per point, today's behavior).	1
seed	int	Random seed.	42

Returns:

Type	Description
dict[str, tuple[NDArray[floating[Any]], NDArray[floating[Any]]]]	Dictionary mapping observable name to an (S1, ST) tuple, each
dict[str, tuple[NDArray[floating[Any]], NDArray[floating[Any]]]]	an array of one value per continuous factor (in the order they
dict[str, tuple[NDArray[floating[Any]], NDArray[floating[Any]]]]	appear in factors).

Raises:

Type	Description
ValueError	If no continuous factors are provided, or n_reps is less than 1.

Source code in src/trade_study/design.py

def sobol_indices(
    run_fn: Callable[[dict[str, Any]], dict[str, float]],
    factors: list[Factor],
    *,
    n_samples: int = 100,
    n_reps: int = 1,
    seed: int = 42,
) -> dict[str, tuple[NDArray[np.floating[Any]], NDArray[np.floating[Any]]]]:
    """Sobol first- and total-order sensitivity indices.

    Like ``screen(method="sobol")``, but returns both S1 (first-order) and
    ST (total-order) per observable instead of discarding ST. ``ST - S1``
    is the standard way to detect interaction effects: a factor with small
    S1 but large ST is interacting with other factors rather than acting
    independently -- information ``screen()`` can't surface at all, since
    it only keeps S1.

    Args:
        run_fn: Callable that takes a config dict and returns a dict of
            observable name -> scalar score.
        factors: Factor list to analyze; only continuous factors are
            varied (as in ``screen()``).
        n_samples: Base sample size *N* for the second-order Saltelli design;
            total configurations are *N* x (2 x num_vars + 2). Physical
            ``run_fn`` calls multiply this by ``n_reps``.
        n_reps: Replicate draws averaged into each sampled point's scores
            before SALib sees them -- see ``screen()``'s ``n_reps``
            for the full rationale and the ``run_fn`` opt-in convention.
            Default 1 (a single draw per point, today's behavior).
        seed: Random seed.

    Returns:
        Dictionary mapping observable name to an ``(S1, ST)`` tuple, each
        an array of one value per continuous factor (in the order they
        appear in ``factors``).

    Raises:
        ValueError: If no continuous factors are provided, or ``n_reps``
            is less than 1.
    """
    from SALib.analyze import sobol as sobol_analyze

    continuous = [f for f in factors if f.factor_type == FactorType.CONTINUOUS]
    if not continuous:
        msg = "Screening requires at least one continuous factor"
        raise ValueError(msg)
    if n_reps < 1:
        msg = f"n_reps must be >= 1; got {n_reps}"
        raise ValueError(msg)

    problem = _salib_problem(continuous)
    results_by_obs = _sobol_sample_and_evaluate(
        run_fn, problem, n_samples, n_reps, seed
    )

    indices: dict[str, tuple[NDArray[np.floating[Any]], NDArray[np.floating[Any]]]] = {}
    for obs_name, vals in results_by_obs.items():
        si = sobol_analyze.analyze(problem, np.array(vals), seed=seed)
        s1 = np.asarray(si["S1"], dtype=np.float64)
        st = np.asarray(si["ST"], dtype=np.float64)
        indices[obs_name] = (s1, st)

    return indices

trade_study.unit_to_value(factor, unit)

Map a unit-interval coordinate to a continuous factor's value.

Parameters:

Name	Type	Description	Default
factor	Factor	Continuous factor with bounds.	required
unit	float	Coordinate in [0, 1].	required

Returns:

Type	Description
float	lo + unit * (hi - lo), or its geometric counterpart for a
float	log-scale factor.

Source code in src/trade_study/design.py

def unit_to_value(factor: Factor, unit: float) -> float:
    """Map a unit-interval coordinate to a continuous factor's value.

    Args:
        factor: Continuous factor with bounds.
        unit: Coordinate in [0, 1].

    Returns:
        ``lo + unit * (hi - lo)``, or its geometric counterpart for a
        log-scale factor.
    """
    lo, hi = _sampling_bounds(factor)
    value = lo + unit * (hi - lo)
    return float(np.exp(value)) if factor.log_scale else float(value)

trade_study.value_to_unit(factor, value)

Map a continuous factor's value to its unit-interval coordinate.

Parameters:

Name	Type	Description	Default
factor	Factor	Continuous factor with bounds.	required
value	float	Factor value within its bounds.	required

Returns:

Type	Description
float	The inverse of :func:unit_to_value.

Source code in src/trade_study/design.py

def value_to_unit(factor: Factor, value: float) -> float:
    """Map a continuous factor's value to its unit-interval coordinate.

    Args:
        factor: Continuous factor with bounds.
        value: Factor value within its bounds.

    Returns:
        The inverse of :func:`unit_to_value`.
    """
    lo, hi = _sampling_bounds(factor)
    coordinate = float(np.log(value)) if factor.log_scale else float(value)
    return (coordinate - lo) / (hi - lo)