Risk Quantification

After solving a scenario program, the complexity k (number of support constraints) and the sample size N are used to compute distribution-free bounds on the probability of constraint violation — the risk — without requiring any knowledge of the underlying probability distribution.

Scenario Approach Background

Consider a decision variable \(x \in \mathcal{X} \subset \mathbb{R}^d\) and an uncertain parameter \(\delta \in \Delta \subset \mathbb{R}^q\) drawn i.i.d. from an unknown probability measure \(\mathbb{P}\). For each scenario \(\delta\), the set \(\mathcal{X}_\delta \subseteq \mathcal{X}\) denotes the region of feasible decisions. The scenario program is:

\[\begin{split}\min_{x \in \mathcal{X}} \quad & c(x) \\ \text{s.t.} \quad & x \in \bigcap_{i=1}^{N} \mathcal{X}_{\delta_i}\end{split}\]

Given the optimal solution \(x_N^*\), the risk quantifies the probability that a newly drawn scenario renders the solution inappropriate:

\[V(x) = \mathbb{P}\!\left\{\delta \in \Delta : x \notin \mathcal{X}_\delta\right\}\]

Since \(\mathbb{P}\) is unknown, the scenario approach provides distribution-free bounds on \(V(x_N^*)\) using only the complexity \(s_N^*\) — the cardinality of the smallest irreducible subset of scenarios that fully determines the solution.

Support Constraints and Complexity

A support list is a sublist \((\delta_{i_1}, \ldots, \delta_{i_k})\) of the scenarios such that:

  1. Solving the program with only these scenarios yields the same solution \(x_N^*\).

  2. The sublist is irreducible: removing any single scenario changes the solution.

The complexity \(s_N^*\) is the minimal cardinality among all support lists. This is the value k computed by get_support() and passed to quantify_risk().

Risk Bounds

Upper Bound (Always Valid)

The upper bound on risk holds under the sole assumption of consistency (which is guaranteed for all convex optimization problems handled by Scen-Opt). For a prescribed confidence level \(1 - \beta\), the risk certificate of Garatti and Campi [GC2025] states:

\[\mathbb{P}^N\!\left\{V(x_N^*) > \epsilon(s_N^*)\right\} \leq \beta\]

where \(\epsilon(k) = 1 - t(k)\) and \(t(k) \in (0,1)\) is the unique solution of

\[\frac{\beta}{N} \sum_{i=k}^{N-1} \binom{i}{k} t^{i-k} \;-\; \binom{N}{k} t^{N-k} \;=\; 0, \qquad k = 0, 1, \ldots, N-1\]

with the convention \(\epsilon(N) = 1\).

Two-Sided Bounds (Requires Non-Degeneracy)

Under the additional assumption of non-degeneracy (there exists a unique support list), a two-sided certificate from Garatti and Campi [GC2022] provides both lower and upper bounds. For \(k = 0, 1, \ldots, N-1\), consider the polynomial equation in \(t\):

\[\binom{N}{k} t^{N-k} \;-\; \frac{\beta}{2N} \sum_{i=k}^{N-1} \binom{i}{k} t^{i-k} \;-\; \frac{\beta}{6N} \sum_{i=N+1}^{4N} \binom{i}{k} t^{i-k} \;=\; 0\]

This admits exactly two solutions in \([0, +\infty)\), denoted \(\underline{t}(k) \leq \overline{t}(k)\). For \(k = N\), a separate equation yields a single solution \(\overline{t}(N)\), with \(\underline{t}(N) = 0\). The bounds are:

\[\underline{\epsilon}(k) = \max\!\big\{0,\; 1 - \overline{t}(k)\big\}, \qquad \overline{\epsilon}(k) = 1 - \underline{t}(k)\]

and the two-sided risk certificate reads:

\[\mathbb{P}^N\!\left\{ \underline{\epsilon}(s_N^*) \;\leq\; V(x_N^*) \;\leq\; \overline{\epsilon}(s_N^*) \right\} \;\geq\; 1 - \beta\]

Scen-Opt computes these bounds numerically via bisection using the regularized incomplete beta function, following the procedure described in [CGC2023].

Warning

The lower bound \(\underline{\epsilon}\) is valid only when the non-degeneracy assumption holds. If Scen-Opt detects degeneracy during active constraint identification (indicated by the degeneracy flag in the solver output), the lower bound should be disregarded and only the upper bound \(\overline{\epsilon}\) used. While Scen-Opt attempts to detect degeneracy automatically, such checks cannot cover out-of-sample scenarios — it is the user’s responsibility to assess whether the lower bound remains applicable.

Tip

A common choice is \(\beta = 10^{-6}\). Smaller values of \(\beta\) give wider bounds but higher confidence.

API Reference

src.Risk.quantify_risk(k, N, beta)[source]

Compute scenario approach risk bounds on constraint violation probability.

Uses the theory of Campi and Garatti with the regularized incomplete beta function to compute distribution-free bounds on the probability of out-of-sample constraint violation via bisection.

The true violation probability satisfies:

\[\varepsilon \;\in\; [\varepsilon_L,\; \varepsilon_U] \quad \text{with confidence at least } 1 - \beta\]
Parameters:

  • k (int) – Cardinality of the support list from the scenario optimization.

  • N (int) – Number of sampled scenarios.

  • beta (float) – Confidence parameter (e.g. 1e-6 for high confidence).

Returns:

  • epsL (float) – Lower bound on the constraint violation probability.

  • epsU (float) – Upper bound on the constraint violation probability.

Raises:

ValueError – If N < 1, k is not between 0 and N, or beta is not strictly between 0 and 1.

Notes

The incomplete beta function betainc used here follows the JAX/SciPy argument convention, which differs from MATLAB (the first argument appears last in JAX).

References

[GC2025]

S. Garatti and M. C. Campi, “Non-convex scenario optimization,” Mathematical Programming, vol. 209, no. 1, pp. 557–608, 2025.

[GC2022]

S. Garatti and M. C. Campi, “Risk and complexity in scenario optimization,” Mathematical Programming, vol. 191, no. 1, pp. 243–279, 2022.

[CGC2023]

M. C. Campi and S. Garatti, “Compression, generalization and learning,” Journal of Machine Learning Research, vol. 24, no. 339, pp. 1–74, 2023.