picos.expressions.uncertain¶
Expression types with an explicit representation of data uncertainty.
Outline¶
Classes¶
A conic description of a |
|
A moment-uncertain description of a random perturbation parameter. |
|
A parameter that can be used to describe (random) noise in data. |
|
Base class for uncertain perturbation sets and distributions. |
|
Base class for random convex or concave piecewise linear expressions. |
|
The maximum over a set of random affine expressions. |
|
The minimum over a set of random affine expressions. |
|
A scenario description of a |
|
A multidimensional uncertain affine expression. |
|
Primary base class for uncertainty affected expression types. |
|
Euclidean or Frobenius norm of an uncertain affine expression. |
|
Squared Euclidean or Frobenius norm of an uncertain affine expression. |
|
Represents perturbation in an Euclidean or Frobenius unit norm ball. |
|
A wasserstein ambiguity set centered at a discrete distribution. |
Exceptions¶
Computing a worst-case (expected) value is hard and not supported. |
Submodules¶
Implements |
|
Implements |
|
Implements |
|
Implements |
|
Implements a parameterization for (random) noise in data. |
|
Implements |
|
Implements |
|
Implements |
|
Implements |
|
Implements the |
Classes¶
ConicPerturbationSet¶
-
class
picos.expressions.uncertain.
ConicPerturbationSet
(parameter_name, shape=(1, 1))[source]¶ Bases:
picos.expressions.uncertain.perturbation.PerturbationUniverse
A conic description of a
Perturbation
.- Definition
An instance
of this class defines a perturbation parameter
where
,
,
,
,
and
is a (product) cone for some
.
- Usage
Obtaining
is done in a number of steps:
Create an instance of this class (see
__init__
).Access
element
to obtain a regular, freshRealVariable
representing.
Define
through any number of regular PICOS constraints that depend only on
and that have a conic representation by passing the constraints to
bound
.Call
compile
to obtain a handle to thePerturbation
.
You can now use
to build instances of
UncertainAffineExpression
and derived constraint types.
It is best practice to assign
to a Python variable and overwrite that variable with
on compilation.
Alternatively, you can obtain a compiled
from the factory method
from_constraints
and accessvia
parameter
.- Example
>>> from picos import Constant, Norm, RealVariable >>> from picos.uncertain import ConicPerturbationSet >>> S = ConicPerturbationSet("P", (4, 4)) >>> P = S.element; P # Obtain a temporary parameter to describe S with. <4×4 Real Variable: P> >>> S.bound(Norm(P, float("inf")) <= 1) # Confine each element to [-1,1]. >>> S.bound(Norm(P, 1) <= 4); S # Allow a perturbation budget of 4. <4×4 Conic Perturbation Set: {P : ‖P‖_max ≤ 1 ∧ ‖P‖_sum ≤ 4}> >>> P = S.compile(); P # Compile the set and obtain the actual parameter. <4×4 Perturbation: P> >>> A = Constant("A", range(16), (4, 4)) >>> U = A + P; U # Define an uncertain data matrix. <4×4 Uncertain Affine Expression: A + P> >>> x = RealVariable("x", 4) >>> U*x # Define an uncertain affine expression. <4×1 Uncertain Affine Expression: (A + P)·x>
-
__init__
(parameter_name, shape=(1, 1))[source]¶ Create a
ConicPerturbationSet
.
-
bound
(constraint)[source]¶ Add a constraint that bounds
.
The constraints do not need to be conic but they need to have a conic representation, which may involve any number of auxiliary variables. For instance, given a constant uncertainty budget
, you may add the bound
(via
picos.Norm(t, 1) <= b
) which can be represented in conic form asThe auxiliary variable
then plays the role of (a slice of)
in the formal definition of
.
When you are done adding bounds, you can obtain
using
compile
.- Raises
RuntimeError – If the set was already compiled.
-
compile
(validate_feasibility=False)[source]¶ Compile the set and return
.
Internally, this computes the matrices
and
, the vector
and the (product) cone
.
- Parameters
validate_feasibility (bool) – Whether to solve the feasibility problem associated with the bounds on
to verify that
is nonempty.
- Returns
An instance of
Perturbation
.- Raises
RuntimeError – If the set was already compiled.
TypeError – If the bound constraints could not be put into conic form.
ValueError – If
could not be verified to be nonempty (needs
validate_feasibility=True
).
-
classmethod
from_constraints
(parameter_name, *constraints)[source]¶ Create a
ConicPerturbationSet
from constraints.The constraints must concern a single regular decision variable that plays the role of the
element
. This variable is not stored or modified and can be reused in a different context.
- Parameters
parameter_name (str) – Name of the parameter that lives in the set.
constraints – A parameter sequence of constraints that concern a single regular decision variable whose internal vectorization is trivial (its dimension must match the product over its shape) and that have a conic representation.
- Raises
ValueError – If the constraints do not all concern the same single variable.
TypeError – If the variable uses a nontrivial vectorization format or if the constraints do not all have a conic representation.
- Example
>>> from picos.expressions.uncertain import ConicPerturbationSet >>> from picos import RealVariable >>> x = RealVariable("x", 4) >>> T = ConicPerturbationSet.from_constraints("t", abs(x) <= 2, x >= 0) >>> print(T) {t : ‖t‖ ≤ 2 ∧ t ≥ 0} >>> print(repr(T.parameter)) <4×1 Perturbation: t>
-
worst_case
(scalar, direction)[source]¶ Implement for
PerturbationUniverse
.
-
property
A
¶ The compiled matrix
.
- Raises
RuntimeError – If the set has not been compiled.
-
property
B
¶ The compiled matrix
or
None
if.
- Raises
RuntimeError – If the set has not been compiled.
-
property
K
¶ The compiled (product) cone
.
- Raises
RuntimeError – If the set has not been compiled.
-
property
c
¶ The compiled vector
.
- Raises
RuntimeError – If the set has not been compiled.
-
property
distributional
¶ Implement for
PerturbationUniverse
.
-
property
element
¶ The perturbation element
describing the set.
This is a regular
RealVariable
that you can use to create constraints to pass tobound
. You can then obtain the “actual” perturbation parameterto use in expressions alongside your decision variaiables using
compile
.Warning
If you use this object instead of
parameter
to define a decision problem then it will act as a regular decision variable, which is probably not what you want.- Raises
RuntimeError – If the set was already compiled.
-
property
ellipsoidal
¶ Whether the perturbation set is an ellipsoid.
If this is true, then a
unit_ball_form
is available.
-
property
parameter
¶ The perturbation parameter
living in the set.
This is the object returned by
compile
.- Raises
RuntimeError – If the set has not been compiled.
-
property
unit_ball_form
¶ A recipe to repose from ellipsoidal to unit norm ball uncertainty.
If the set is
ellipsoidal
, then this is a pair(U, M)
whereU
is aUnitBallPerturbationSet
andM
is a dictionary mapping the oldparameter
to an affine expression of the new parameter that can represent the old parameter in an expression (seereplace_mutables
). The mappingM
is empty if and only if the perturbation set is already an instance ofUnitBallPerturbationSet
.If the uncertainty set is not ellipsoidal, then this is
None
.See also
SOCConstraint.unit_ball_form
.- Example
>>> from picos import Problem, RealVariable, sum >>> from picos.uncertain import ConicPerturbationSet >>> # Create a conic perturbation set and a refinement recipe. >>> T = ConicPerturbationSet("t", (2, 2)) >>> T.bound(abs(([[1, 2], [3, 4]] ^ T.element) + 1) <= 10) >>> t = T.compile() >>> U, mapping = T.unit_ball_form >>> print(U) {t' : ‖t'‖ ≤ 1} >>> print(mapping) {<2×2 Perturbation: t>: <2×2 Uncertain Affine Expression: t(t')>} >>> # Define and solve a conically uncertain LP. >>> X = RealVariable("X", (2, 2)) >>> P = Problem() >>> P.set_objective("max", sum(X)) >>> _ = P.add_constraint(X + 2*t <= 10) >>> print(repr(P.parameters["t"].universe)) <2×2 Conic Perturbation Set: {t : ‖[2×2]⊙t + [1]‖ ≤ 10}> >>> _ = P.solve(solver="cvxopt") >>> print(X) [-8.00e+00 1.00e+00] [ 4.00e+00 5.50e+00] >>> # Refine the problem to a unit ball uncertain LP. >>> Q = Problem() >>> Q.set_objective("max", sum(X)) >>> _ = Q.add_constraint(X + 2*mapping[t] <= 10) >>> print(repr(Q.parameters["t'"].universe)) <2×2 Unit Ball Perturbation Set: {t' : ‖t'‖ ≤ 1}> >>> _ = Q.solve(solver="cvxopt") >>> print(X) [-8.00e+00 1.00e+00] [ 4.00e+00 5.50e+00]
MomentAmbiguitySet¶
-
class
picos.expressions.uncertain.
MomentAmbiguitySet
(parameter_name, shape, nominal_mean=0, nominal_covariance='I', alpha=0, beta=1, sample_space=None)[source]¶ Bases:
picos.expressions.uncertain.perturbation.PerturbationUniverse
A moment-uncertain description of a random perturbation parameter.
- Model of uncertainty
As a distributional ambiguity set, an instance
of this class
represents a safety region for a partially known (ambiguous) probability distribution
and
provides a random, ambiguously distributed perturbation parameter
that can be used to define worst-case-expectation expressions of the form
for a selection of functions
and a decision variable
.
- Definition
Formally, this class can describe ambiguity sets of the form
where
is the set of all Borel probability measures on
,
the sample space
bounds the support of
and may be given as either
or as an
-dimensional ellipsoid,
and
with
are the nominal mean and covariance matrix of
, respectively, and
and
are meta-parameters bounding the uncertainty concerning the mean and the covariance matrix, respectively.
Unless
, this class can also represent the limit cases of
and
denoting an unbounded mean and covariance matrix, respectively.
Note
While
denotes that the nominal mean
is certain, there is a subtle difference between setting
on the one hand and assuming a certain form for the covariance matrix on the other hand: In the former case, the worst case covariance matrix may be Lowener smaller than the nominal one. Setting a lower bound on the covarianve matrix is computationally difficult and not supported.
- Supported functions
A squared norm
where
is biaffine in
and
. This can be written as
abs(A)**2
in Python.A convex piecewise linear function
where
is biaffine in
and
for all
. This can be written as
picos.max([a_1, ..., a_k])
in Python.A concave piecewise linear function
where
is biaffine in
and
for all
. This can be written as
picos.min([a_1, ..., a_k])
in Python.
- Example
We show that for unbounded mean and support (i.e.
and
) and for a finite samples space
, this distributional ambiguity set can be used in the context of classical (non-distributional) robust optimization applied to least squares problems.
>>> from picos import Constant, diag, Ellipsoid, Problem, RealVariable >>> from picos.uncertain import ConicPerturbationSet, MomentAmbiguitySet >>> import numpy >>> numpy.random.seed(1) >>> # Setup data and variables of the nominal least squares problem. >>> n = 3 >>> A = Constant("A", numpy.random.random((n, n))) >>> b = Constant("b", numpy.random.random(n)) >>> x = RealVariable("x", n) >>> # Setup an ellipsoid S bounding the uncertainty in both models. >>> N = n**2 >>> S = Ellipsoid(N, diag(range(1, N + 1)), range(N)) >>> # Define and solve both robust models. >>> U1 = ConicPerturbationSet.from_constraints( ... "Y", RealVariable("Y", N) << S) >>> U2 = MomentAmbiguitySet("Z", N, alpha=None, beta=None, sample_space=S) >>> Y = U1.parameter.reshaped((n, n)) >>> Z = U2.parameter.reshaped((n, n)) >>> P1, P2 = Problem(), Problem() >>> P1.objective = "min", abs((A + Y)*x - b) >>> P2.objective = "min", abs((A + Z)*x - b)**2 >>> _ = P1.solve(solver="cvxopt") >>> x1 = ~x # Save current value of x as a constant PICOS expression x1. >>> _ = P2.solve(solver="cvxopt") >>> x2 = ~x >>> # Verify that both models yield the same robust regression vector. >>> x1.equals(x2, relTol=1e-4) True
-
__init__
(parameter_name, shape, nominal_mean=0, nominal_covariance='I', alpha=0, beta=1, sample_space=None)[source]¶ Create a
MomentAmbiguitySet
.- Parameters
parameter_name (str) – Name of the random parameter
.
shape (anything recognized by
load_shape
) – Shape of. Must characterize a column vector. If
None
, then the shape is inferred from the nominal mean.nominal_mean (anything recognized by
load_data
) – The nominal meanof the ambiguous distribution
.
nominal_covariance (anything recognized by
load_data
) – The nominal covariance matrixof the ambiguous distribution
.
alpha (float) – The parameter
bounding the uncertain mean. The values
None
andfloat("inf")
denote an unbounded mean.beta (float) – The parameter
bounding the uncertain covariance matrix. The values
None
andfloat("inf")
denote unbounded covariances.sample_space (
None
orTheField
orEllipsoid
) – The sample space. If this is
None
or an instance ofTheField
(i.e.), then the support of
is unbounded. If this is an
-dimensional instance of
Ellipsoid
, then the support ofis a subset of that ellipsoid.
-
property
dim
¶ The dimension
of the sample space.
-
property
distributional
¶ Implement for
PerturbationUniverse
.
-
property
nominal_covariance
¶ The nominal covariance matrix
.
-
property
nominal_mean
¶ The nominal mean
of the ambiguous distribution.
-
property
parameter
¶ The random perturbation parameter
.
Perturbation¶
-
class
picos.expressions.uncertain.
Perturbation
(universe, name, shape)[source]¶ Bases:
picos.expressions.mutable.Mutable
,picos.expressions.uncertain.uexp_affine.UncertainAffineExpression
A parameter that can be used to describe (random) noise in data.
This is the initial building block for an
UncertainAffineExpression
. In particular, an affine transformation of this parameter represents uncertain data.-
__init__
(universe, name, shape)[source]¶ Create a
Perturbation
.- Parameters
universe (PerturbationUniverse) – Either the set that the perturbation parameter lives in or the distribution according to which the perturbation is distributed.
name (str) – Symbolic string description of the perturbation, similar to a variable’s name.
shape (int or tuple or list) – Algebraic shape of the perturbation parameter.
This constructor is meant for internal use. As a user, you will want to first define a universe (e.g.
ConicPerturbationSet
) for the parameter and obtain the parameter from it.
-
property
universe
¶ The uncertainty universe that the parameter belongs to.
-
PerturbationUniverse¶
-
class
picos.expressions.uncertain.
PerturbationUniverse
[source]¶ Bases:
abc.ABC
Base class for uncertain perturbation sets and distributions.
See
distributional
for a distinction between perturbation sets, random distributions and distributional ambiguity sets, all three of which can be represented by this class.The naming scheme for implementing classes is as follows:
Perturbation sets (robust optimization) end in
PerturbationSet
,random distributions (stochastic programming) end in
Distribution
,distributional ambiguity sets (DRO) end in
AmbiguitySet
.
-
classmethod
make_type
(*args, **kwargs)[source]¶ Create a detailed universe type from subtype parameters.
-
worst_case
(scalar, direction)[source]¶ Find a worst-case realization of the uncertainty for an expression.
- Parameters
scalar (UncertainExpression) – A scalar uncertain expression that depends only on the perturbation
parameter
.direction (str) – Either
"min"
or"max"
, denoting the worst-case direction.
- Returns
A pair where the first element is the worst-case (expeceted) value as a
float
and where the second element is a realization of the perturbation parameter that attains this worst case as afloat
or CVXOPT matrix (orNone
for stochastic uncertainty).- Raises
TypeError – When the function is not scalar.
ValueError – When the function depends on other mutables than exactly the
parameter
.picos.uncertain.IntractableWorstCase – When computing the worst-case (expected) value is not supported, in particular when it would require solving a nonconvex problem.
RuntimeError – When the computation is supported but fails.
-
abstract property
distributional
¶ Whether this is a distribution or distributional ambiguity set.
If this is
True
, then this represents a random distribution (stochastic programming) or an ambiguity set of random distributions (distributionally robust optimization) and any expression that depends on its randomparameter
, when used in a constraint or as an objective function, is understood as a (worst-case) expected value.If this is
False
, then this represents a perturbation set (robust optimization) and any expression that depends on its perturbationparameter
, when used in a constraint or as an objective function, is understood as a worst-case value.
-
abstract property
parameter
¶ The perturbation parameter.
-
property
subtype
¶
-
property
type
¶ Detailed type of a perturbation parameter universe.
RandomExtremumAffine¶
-
class
picos.expressions.uncertain.
RandomExtremumAffine
(expressions)[source]¶ Bases:
picos.expressions.exp_extremum.ExtremumBase
,picos.expressions.uncertain.uexpression.UncertainExpression
,picos.expressions.expression.Expression
Base class for random convex or concave piecewise linear expressions.
Note
Unlike other uncertain expression types, this class is limited to uncertainty of stochastic nature, where using the expression in a constraint or as an objective function implicitly takes the (worst-case) expectation of the expression. Non-stochastic uncertainty is handled within
MaximumConvex
andMinimumConcave
as their behavior, although designed for certain expression types, already encodes the worst-case approach of the robust optimization paradigm.-
__init__
(expressions)[source]¶ Construct a
RandomExtremumAffine
.- Parameters
expressions – A collection of uncertain affine expressions whose uncertainty is of stochastic nature.
-
property
expressions
¶ The expressions under the extremum.
-
property
perturbation
¶ Fast override for
UncertainExpression
.
-
RandomMaximumAffine¶
-
class
picos.expressions.uncertain.
RandomMaximumAffine
(expressions)[source]¶ Bases:
picos.expressions.exp_extremum.MaximumBase
,picos.expressions.uncertain.uexp_rand_pwl.RandomExtremumAffine
The maximum over a set of random affine expressions.
RandomMinimumAffine¶
-
class
picos.expressions.uncertain.
RandomMinimumAffine
(expressions)[source]¶ Bases:
picos.expressions.exp_extremum.MinimumBase
,picos.expressions.uncertain.uexp_rand_pwl.RandomExtremumAffine
The minimum over a set of random affine expressions.
ScenarioPerturbationSet¶
-
class
picos.expressions.uncertain.
ScenarioPerturbationSet
(parameter_name, scenarios, compute_hull=None)[source]¶ Bases:
picos.expressions.uncertain.perturbation.PerturbationUniverse
A scenario description of a
Perturbation
.- Definition
An instance
of this class defines a perturbation parameter
where
is a finite set of scenarios and
denotes the convex hull.
Usually, the scenarios are observed or projected realizations of the uncertain data and
is used to represent the data directly.
- Example
>>> from picos.uncertain import ScenarioPerturbationSet >>> scenarios = [[1, -1], [1, 1], [-1, -1], [-1, 1], [0, 0]] >>> S = ScenarioPerturbationSet("s", scenarios, False); S <2×1 Scenario Perturbation Set: conv({5 2×1 scenarios})> >>> s = S.parameter; s <2×1 Perturbation: s> >>> # Compute largest sum of entries over all points in S. >>> value, realization = S.worst_case(s[0] + s[1], "max") >>> round(value, 4) 2.0 >>> print(realization) [ 1.00e+00] [ 1.00e+00]
-
__init__
(parameter_name, scenarios, compute_hull=None)[source]¶ Create a
ScenarioPerturbationSet
.- Parameters
parameter_name (str) – Name of the parameter that lives in the set.
scenarios (anything recognized by
picos.Samples
) – A collection of data points of same shape representing.
compute_hull (bool) – Whether to use SciPy to compute the convex hull of the data points and discard points in the interior. This can speed up the solution process significantly, in particular when the scenarios come from observations and when the data is low-dimensional. On the other hand, when the given scenarios are known to be on the boundary of their convex hull, then disabling this speeds up initialization of the perturbation set. The default value of
None
meansTrue
when SciPy is available andFalse
otherwise.
-
worst_case
(scalar, direction)[source]¶ Implement for
PerturbationUniverse
.
-
property
distributional
¶ Implement for
PerturbationUniverse
.
-
property
parameter
¶ Implement for
PerturbationUniverse
.
UncertainAffineExpression¶
-
class
picos.expressions.uncertain.
UncertainAffineExpression
(string, shape=(1, 1), coefficients={})[source]¶ Bases:
picos.expressions.uncertain.uexpression.UncertainExpression
,picos.expressions.exp_biaffine.BiaffineExpression
A multidimensional uncertain affine expression.
This expression has the form
where
,
,
,
and
are defined as for the
BiaffineExpression
base class andis an uncertain perturbation parameter confined to (distributed according to) a perturbation set (distribution)
.
If no coefficient matrices defining
and
are provided, then this expression represents uncertain data confined to an uncertainty set
(distributed according to
) where
can be understood as a nominal data value while
quantifies the uncertainty on the data.
-
__init__
(string, shape=(1, 1), coefficients={})[source]¶ Construct an
UncertainAffineExpression
.Extends
exp_biaffine.BiaffineExpression.__init__
.This constructor is meant for internal use. As a user, you will want to first define a universe (e.g.
ConicPerturbationSet
) for aperturbation parameter
and use that parameter as a building block to create more complex uncertain expressions.
-
UncertainExpression¶
-
class
picos.expressions.uncertain.
UncertainExpression
[source]¶ Bases:
object
Primary base class for uncertainty affected expression types.
The secondary base class must be
Expression
or a subclass thereof.Uncertain expressions have a distinct behavior when used to form a constraint or when posed as an objective function. The exact behavior depends on the type of uncertainty involved. If the perturbation parameter that describes the uncertainty is confied to a perturbation set, then the worst-case realization of the parameter is assumed when determining feasibility and optimality. If the perturbation parameter is a random variable (whose distribution may itself be ambiguous), then the constraint or objective implicitly considers the expected value of the uncertain expression (under the worst-case distribution). Uncertain expressions are thus used in the contexts of robust optimization, stochastic programming and distributionally robust optimization.
-
worst_case
(direction)[source]¶ Find a worst-case realization of the uncertainty for the expression.
Expressions that are affected by uncertainty are only partially valued once an optimization solution has been applied. While their decision values are populated with a robust optimal solution, the parameter that controls the uncertainty is not valued unless the user assigned it a particular realization by hand. This method computes a worst-case (expected) value of the expression and returns it together with a realization of the perturbation parameter for which the worst case is attained (or
None
in the case of stochastic uncertainty).For multidimensional expressions, this method computes the entrywise worst case and returns an attaining realization for each entry.
- Parameters
direction (str) – Either
"min"
or"max"
, denoting the worst-case direction.- Returns
A pair
(value, realization)
. For a scalar expression,value
is its worst-case (expected) value as afloat
andrealization
is a realization of theperturbation
parameter that attains this worst case as afloat
or CVXOPT matrix. For a multidimensional expression,value
is a CVXOPT dense matrix denoting the entrywise worst-case values andrealization
is atuple
of attaining realizations corresponding to the expression vectorized in in column-major order. Lastly,realization
isNone
if the expression iscertain
or when its uncertainty is of stochastic nature.- Raises
picos.NotValued – When the decision variables that occur in the expression are not fully valued.
picos.uncertain.IntractableWorstCase – When computing the worst-case (expected) value is not supported, in particular when it would require solving a nonconvex problem.
RuntimeError – When the computation is supported but fails.
-
worst_case_string
(direction)[source]¶ A string describing the expression within a worst-case context.
- Parameters
direction (str) – Either
"min"
or"max"
, denoting the worst-case direction.
-
worst_case_value
(direction)[source]¶ A shorthand for the first value returned by
worst_case
.
-
property
certain
¶ Whether the uncertain expression is actually certain.
-
property
random
¶ Whether the uncertainty is of stochastic nature.
See also
distributional
.
-
property
uncertain
¶ Whether the uncertain expression is in fact uncertain.
-
property
universe
¶ Universe that the perturbation parameter lives in, or
None
.If this is not
None
, then this is the same asperturbation
.:attr:~.perturbation.Perturbation.universe.
-
UncertainNorm¶
-
class
picos.expressions.uncertain.
UncertainNorm
(x)[source]¶ Bases:
picos.expressions.uncertain.uexpression.UncertainExpression
,picos.expressions.expression.Expression
Euclidean or Frobenius norm of an uncertain affine expression.
-
__init__
(x)[source]¶ Construct an
UncertainNorm
.- Parameters
x (UncertainAffineExpression) – The uncertain affine expression to denote the norm of.
-
property
x
¶ Uncertain affine expression under the norm.
-
UncertainSquaredNorm¶
-
class
picos.expressions.uncertain.
UncertainSquaredNorm
(x)[source]¶ Bases:
picos.expressions.uncertain.uexpression.UncertainExpression
,picos.expressions.expression.Expression
Squared Euclidean or Frobenius norm of an uncertain affine expression.
-
__init__
(x)[source]¶ Construct an
UncertainSquaredNorm
.- Parameters
x (UncertainAffineExpression) – The uncertain affine expression to denote the squared norm of.
-
property
x
¶ Uncertain affine expression under the squared norm.
-
UnitBallPerturbationSet¶
-
class
picos.expressions.uncertain.
UnitBallPerturbationSet
(parameter_name, shape=(1, 1))[source]¶ Bases:
picos.expressions.uncertain.pert_conic.ConicPerturbationSet
Represents perturbation in an Euclidean or Frobenius unit norm ball.
This is a
ConicPerturbationSet
with fixed formAfter initialization, you can obtain the parameter using
parameter
.-
property
unit_ball_form
¶ Overwrite
ConicPerturbationSet.unit_ball_form
.
-
property
WassersteinAmbiguitySet¶
-
class
picos.expressions.uncertain.
WassersteinAmbiguitySet
(parameter_name, p, eps, samples, weights=1)[source]¶ Bases:
picos.expressions.uncertain.perturbation.PerturbationUniverse
A wasserstein ambiguity set centered at a discrete distribution.
- Model of uncertainty
As a distributional ambiguity set, an instance
of this class
represents a safety region for a partially known (ambiguous) probability distribution
and
provides a random, ambiguously distributed perturbation parameter
that can be used to define worst-case-expectation expressions of the form
for a selection of functions
and a decision variable
.
- Definition
Formally, this class can describe discrepancy-based ambiguity sets of the form
where discrepancy from the discrete nominal distribution
is measured with respect to the Wasserstein distance of order
,
where
is the set of all Borel probability measures on
for some
,
denotes the set of all couplings of
and
,
for all
are the
samples comprising the support of
,
are weights denoting the nominal probabilitiy mass at
for all
,
denotes the Dirac delta function with unit mass at
for all
and where
controls the radius of the ambiguity set.
- Supported functions
For
:
A convex piecewise linear function
where
is biaffine in
and
for all
. This can be written as
picos.max([a_1, ..., a_k])
in Python.A concave piecewise linear function
where
is biaffine in
and
for all
. This can be written as
picos.min([a_1, ..., a_k])
in Python.
For
:
A squared norm
where
is biaffine in
and
. This can be written as
abs(A)**2
in Python.
-
__init__
(parameter_name, p, eps, samples, weights=1)[source]¶ Create a
WassersteinAmbiguitySet
.- Parameters
parameter_name (str) – Name of the random parameter
.
p (float) – The Wasserstein type/order parameter
.
eps (float) – The Wasserstein ball radius
.
samples (aynthing recognized by
Samples
) – The support of the discrete distributiongiven as the samples
. The original shape of the samples determines the shape of
.
weights – A vector denoting the nonnegative weight (e.g. frequency or probability) of each sample. Its length must match the number of samples provided. The argument will be normalized such that its entries sum to one. Entries of zero will be dropped alongside their associated sample. The default value of
1
denotes the empirical distribution on the samples.
Warning
Duplicate samples are not detected and can impact performance. If duplicate samples are likely, make sure to detect them and encode their frequency in the weight vector.
-
property
distributional
¶ Implement for
PerturbationUniverse
.
-
property
eps
¶ The Wasserstein ball radius
.
-
property
p
¶ The Wasserstein order
.
-
property
parameter
¶ The random perturbation parameter
.
-
property
weights
¶ The sample weights a constant PICOS vector.
Exceptions¶
IntractableWorstCase¶
-
exception
picos.expressions.uncertain.
IntractableWorstCase
[source]¶ Bases:
RuntimeError
Computing a worst-case (expected) value is hard and not supported.
Raised by
worst_case
and methods that depend on it.-
__init__
(*args, **kwargs)¶ Initialize self. See help(type(self)) for accurate signature.
-
__new__
(**kwargs)¶ Create and return a new object. See help(type) for accurate signature.
-