Module index

pythia.basis

Assemble sparse univariate and multivariate basis polynomials.

pythia.basis.multivariate_basis(univariate_bases, indices, partial=None)[source]

Assemble multivariate polynomial basis.

Set the (partial derivative of the) multivariate (product) polynomial basis functions.

Parameters

univariate_bases (list of list of callable) – Univariate basis functions for parameters. Is called by univariate_bases[paramIdx][deg]().
indices (array_like) – Array of multiindices for multivariate basis functions.
partial (list of int) – Number of partial derivatives for each dimension. Length is same as univariate_bases.

Return type

list[Callable]

Returns

List of multivariate product polynomials with univariate degrees as specified in indices.

pythia.basis.normalize_polynomial(weight, basis, param)[source]

Normalize orthogonal polynomials.

Normalize a polynomial of an orthogonal system with respect to the scalar product

\[a(u,v)_\mathrm{pdf} = \int u(p) v(p) \mathrm{pdf}(p) \mathrm{d}p.\]

The normalized polynomial \(\phi_j\) for any given polynomial \(P_j\) is given by \(\phi_j = P_j / \sqrt{c_j}\) for the constant \(c_j = \int \mathrm{pdf}(p) * P_j(p)^2 \mathrm{d}p\).

Parameters

weight (callable) – Probability density function.
basis (list of numpy.polynomial.Polynomial) – Polynomials to normalize w.r.t. weight.
param (pythia.parameter.Parameter) – Parameter used for distribution and domain information.

Return type

list[Polynomial]

Returns

List of normalized univariate polynomials.

pythia.basis.set_hermite_basis(param, deg)[source]

Generate list of probabilists Hermite polynomials.

Generate the Hermite Polynomials up to certain degree according to the mean and variance of the specified parameter.

Parameters

param (pythia.parameters.Parameter) – Parameter for basis function. Needs to be normal distributed.
deg (int) – Maximum degree for polynomials.

Return type

list[Polynomial]

Returns

List of probabilists Hermite polynomials up to (including) degree specified in deg.

pythia.basis.set_jacobi_basis(param, deg)[source]

Generate list of Jacobi polynomials.

Generate the Jacobi Polynomials up to certain degree on the interval and DoFs specified by the parameter.

Note

The Jacobi polynomials have leading coefficient 1.

Parameters

param (pythia.parameters.Parameter) – Parameter for basis function. Needs to be Beta-distributed.
deg (int) – Maximum degree for polynomials.

Return type

list[Polynomial]

Returns

List of Jacobi polynomials up to (including) degree specified in deg.

pythia.basis.set_laguerre_basis(param, deg)[source]

Generate list of Leguerre polynomials.

Generate the generalized Laguerre polynomials up to certain degree on the interval and DoFs specified by the parameter.

Parameters

param (pythia.parameters.Parameter) – Parameter for basis function. Needs to be Gamma-distributed.
deg (int) – Maximum degree for polynomials.

Return type

list[Polynomial]

Returns

List of Laguerre polynomials up to (including) degree specified in deg.

pythia.basis.set_legendre_basis(param, deg)[source]

Generate list of the Legendre Polynomials.

Generate the Legendre Polynomials up to certain degree on the interval specified by the parameter.

Parameters

param (pythia.parameters.Parameter) – Parameter for basis function. Needs to be uniformly distributed.
deg (int) – Maximum degree for polynomials.

Return type

list[Polynomial]

Returns

List of Legendre polynomials up to (including) degree specified in deg.

pythia.basis.transform_interval(points, source_interval, target_interval)[source]

Affinely transform points from source_inverval to target_interval.

Parameters

points (array_like) – Points lying in the source_interval.
source_interval (array_like) – Source interval.
target_interval (array_like) – Target interval.

Return type

ndarray

Returns

Transformed points lying in the target_interval.

pythia.basis.univariate_basis(params, degs)[source]

Assemble a univariate polynomial basis.

Set polynomial basis up to deg for each parameter in params according to the parameter distribution and area of definition.

Parameters

params (list of pythia.parameter.Parameter) – Parameters to compute univariate basis function for.
degs (array_like) – Max. degrees of univariate polynomials for each parameter.

Return type

list[list[Polynomial]]

Returns

List of normalized univariate polynomials w.r.t. parameter domain and distribution up to specified degree for each parameter in params.

pythia.chaos

Sample-based computation of polynomial chaos expansion.

class pythia.chaos.PolynomialChaos(params, index_set, x_train, w_train, y_train, coefficients=None)[source]

Bases: object

Computation of sparse polynomial chaos expansion.

Parameters

params (list of pt.parameter.Parameter) – List of stochastic parameters.
index_set (pt.index.IndexSet) – Index set for sparse polynomial chaos expansion.
x_train (array_like) – Parameter realizations for training.
weights (array_like) – Regression weights for training.
fEval (array_like) – Function evaluation for training.
coefficients (array_like, default=None) – Polynomial expansion coefficients. If given, the coefficients are not computed during initiation. This can be used to load a chaos expansion.

eval(x, partial=None)[source]

Evaluate the (partial derivative of the) PC approximation.

Parameters

x (np.ndarray) – Parameter realizations in which the approximation is evaluated.
partial (list[int] | dict | None, optional) – Number of derivatives for each parameter component. If a list is given, length has to be the number of parameters. Ordering of list is according to self.parameters. If a dict is given, keys have to be subset of parameter names.

Return type

ndarray

Returns

Evaluation of polynomial expansion in x values.

Examples

Given two parameters \(x_1\) and \(x_2\)

>>> param1 = pt.parameter.Parameter("x1", [-1, 1], "uniform")
>>> param2 = pt.parameter.Parameter("x2", [-1, 1], "uniform")

a corresponding polynomial chaos approximation for a function \(f: (x_1,x_1) \mapsto y\)

>>> surrogate = pt.chaos.PolynomialChaos([param1, param2], ...)

and an array the surrogate should be evaluated in

>>> x_test = np.random.uniform(-1, 1, (1000, 2))

we can evaluate the surrogate with

>>> y_approx = surrogate.eval(x_test)

To obtain partial a partial derivative of the approximation, e.g., \(\frac{\partial^2f}{\partial x_2^2}\), specify a list

>>> y_approx = surrogate.eval(x_test, partial=[0, 2])

or a dictionary with parameter names and number of partial derivates

>>> y_approx = surrogate.eval(x_test, partial={'x2':2})

property mean: ndarray: Mean of the PC expansion.

property std: ndarray: Standard deviation of the PC expansion.

property var: ndarray: Variance of the PC expansion.

pythia.chaos.assemble_indices(enum_idx, sobol_tuples, max_terms)[source]

Compute automatic choice of multiindices.

Parameters

enum_idx (np.ndarray) – Sorted enumeration indices according to magnitude of Sobol indices.
sobol_tuples (list of tuple) – List of Sobol subscript tuples.
max_terms (int) – Maximum number of expansion terms.

Returns

indices – Array of (sparse) optimal multiindices.

Return type

np.ndarray

pythia.chaos.doerfler_marking(values, idx=None, threshold=0.9)[source]

Dörfler marking for arbitrary values.

Parameters

values (array_like) – Values for the Dörfler marking.
idx (list of int, optional) – List of indices associated with the entries of values. If None, this is set to range(len(values)).
threshold (float, default=0.9) – Threshold paramter for Dörfler marking.

Return type

tuple[ndarray, ndarray, int]

Returns

idx_reordered – Reordered indices given by idx. Ordered from largest to smallest value.
ordered_values – Reordered values. Ordered from largest to smallest.
marker – Threshold marker such that sum(values[:marker]) > threshold * sum(values).

pythia.chaos.find_optimal_indices(params, x_train, w_train, y_train, max_terms=0, threshold=0.001)[source]

Compute optimal multiindices of polynomial chaos expansion.

Heuristical approach to compute almost optimal multiindices for a polynomial chaos expansion based on an estimate of the Sobol index values.

Parameters

params (list of pythia.Parameters.Parameter) – Random parameters of the problem.
x_train (array_like) – Sample points for training
w_train (array_like) – Weights for training.
y_train (array_like) – Function evaluations for training.
max_terms (int, default=0) – Maximum number of expansion terms. Number of expansion terms is chosen automatically for max_terms=0.
threshold (float, default=1e-03) – Truncation threshold for Sobol indices. Smallest Sobol values with sum less then threshold are ignored.

Return type

tuple[ndarray, ndarray]

Returns

indices – Array with multiindices.
sobol – Crude intermediate approximation of Sobol indices.

Notes

To find reasonable candidates for the sparse polynomial chaos expansion, first an expansion with a large simplex index set is computed. The simplex index set uses the same maximum dimension in each component and is designed to have at least max_terms many elements. With this index set a polynomial chaos expansion is computed. The computed Sobol indices are then ordered and the largest contributions are collected by a Dörfler marking strategy. Then a new index set is assembled by including a downward closed subset of polynomial chaos coefficient indices for each selected Sobol index tuple. The number of chosen indices for each selected Sobol index tuple is weighted by the respective Sobol index value.

pythia.chaos.get_gram_batchsize(dim, save_memory=538445312.5)[source]

Compute memory allocation batch sizes for information matrix.

Compute the maximal number of samples in each batch when assembling the information matrix to be maximally memory efficient and avoid OutOfMemory errors.

Parameters

dim (int) – Number of rows/columns of information matrix.
save_memory (int, default=3*1025/2) – Memory (in bytes), that should be kept free. The default is equivalent to 512 MB.

Return type

int

Returns

Batchsize for assembling of information matrix.

pythia.index

Create, manipulate and store information about multiindices.

class pythia.index.IndexSet(indices)[source]

Bases: object

Generate index set object for sparse PC expansion.

A general polynomial chaos expansion of a function \(f\colon\Gamma\subset\mathbb{R}^M\to\mathbb{R}^J\) with \(y\sim\pi\) is given by

\[f(y) = \sum_{\mu\in\mathbb{N}_0^M} \mathbf{f}[\mu]P_{\mu}(y) \quad\mbox{for}\quad \mathbf{f}[\mu] = \int_\Gamma f(y)P_\mu(y)\ \mathrm{d}y,\]

where \(\mu\) is a multiindex, \(\mathbf{f}[\mu]\in\mathbb{R}^J\) is a coefficient vector and \(\{P_\mu\}_{\mu\in\mathbb{N}_0^M}\) is an orthonormal basis in \(L^2(\Gamma,\pi)\). To approximate the infinite expansion choose an index set \(\Lambda\subset\mathbb{N}_0^M\) of multiindices and consider

\[f(y) \approx \sum_{\mu\in\Lambda} \mathbf{f}[\mu]P_{\mu}(y),\]

Parameters: indices (np.ndarray) – Array of multiindices with shape (#indices, param dim).

Examples

Create the sparse index set

\[\Lambda = \{ (0,0), (1,0), (2,0), (0,1) \} \subset \mathbb{N}_0^2\]

>>> import pythia as pt
>>> indices = np.array([[0, 0], [1, 0], [2, 0], [0, 1]], dtype=int)
>>> index_set = pt.index.IndexSet(indices)

get_index_number(indices)[source]

Get enumeration number of indices.

Get the row indices of the given multiindices such that self.indices[rows] = indices.

Parameters: indices (np.ndarray) – Indices to get the number of.
Return type: ndarray
Returns: Array containing the enumeration numbers of the indices.

get_sobol_tuple_number(sobol_tuples)[source]

Get enumeration indices of Sobol tuples.

Parameters: sobol_tuples (list of tuple) – List of Sobol tuples.
Return type: ndarray
Returns: Array containing the enumeration number of the Sobol tuples.

index_to_sobol_tuple(indices)[source]

Map array of indices to their respective Sobol tuples.

Parameters: indices (np.ndarray) – Array of multiindices.
Return type: list[tuple]
Returns: List of Sobol tuples.

sobol_tuple_to_indices(sobol_tuples)[source]

Map Sobol tuples to their respective indices.

Parameters: sobol_tuples (tuple or list of tuple) – List of Sobol tuples.
Return type: list[ndarray]
Returns: List of index arrays for each given Sobol tuple.

pythia.index.intersection(index_list)[source]

Intersect list of multiindex sets.

Given sparse index sets \(\Lambda_1, \dots, \Lambda_N\), compute \(\Lambda=\Lambda_1\cap\dots\cap\Lambda_N\).

Parameters: index_list (list[np.ndarray]) – List of index sets.
Return type: ndarray
Returns: Intersection of index sets.

pythia.index.lq_bound_set(dimensions, bound, q=1.0)[source]

Create set of multiindices with bounded \(\ell^q\)-norm.

For given dimensions \(d \in \mathbb{N}^M\), bound \(b \in \mathbb{R}_{>0}\) and norm factor \(q \in \mathbb{R}_{>0}\), the \(\ell^q\)-norm index set is given by

\[\Lambda = \{ \mu\in [d_1]\times\dots\times [d_M] \ \vert\ \Vert \mu \Vert_{\ell^q} \leq b\},\]

where \([d_m]=\{0, \dots, d_m-1\}\) and

\[\Vert \mu \Vert_{\ell^q} = \Bigl(\sum_{m=1}^M \mu_m^q\Bigr)^{\frac{1}{q}}.\]

Parameters

dimensions (Sequence[int]) – Dimensions for each component, i.e., indices from 0 to dimension-1.
bound (float) – Bound for the \(\ell^q\)-norm.
q (float, optional) – Norm factor.

Return type

ndarray

Returns

Array with all possible multiindices with bounded \(\ell^q\)-norm.

Examples

>>> pt.index.lq_bound_set([5, 5], 4, 0.5)
array([[0, 0],
       [0, 1],
       [1, 0],
       [0, 2],
       [1, 1],
       [2, 0],
       [0, 3],
       [3, 0],
       [0, 4],
       [4, 0]])

pythia.index.set_difference(indices, subtract)[source]

Set difference of two index arrays.

Given two sparse index sets \(\Lambda_1\) and \(\Lambda_2\), compute \(\Lambda=\Lambda_1\setminus\Lambda_2\).

Parameters

indices (np.ndarray) – Index array multiindices are taken out of.
subtract (np.ndarray) – Indices that are taken out of the original set.

Return type

ndarray

Returns

Set difference of both index arrays.

pythia.index.simplex_set(dimension, maximum)[source]

Create a simplex index set.

For given dimension \(M\in\mathbb{N}\) and maximum \(d\in\mathbb{N}\) the simplex index set is given by

\[\Lambda = \{ \mu\in\mathbb{N}_0^M \ \vert\ \sum_{m=1}^M \mu_m \leq d\}.\]

Notes

Limiting the absolute value of the multiindices creates a simplex in \(\mathbb{N}_0^M\), which motivates the name of the function. As an example, in two dimensions this gives us points inside a triangle limited by the axes and the line \(x_1 + x_2 = d\).

Parameters

dimension (int) – Dimension of the multiindices.
maximum (int) – Maximal sum value for the multiindices.

Return type

ndarray

Returns

Array with all possible multiindices in simplex set.

Examples

>>> pt.index.simplex(2, 2)
array([[0, 0],
       [0, 1],
       [1, 0],
       [0, 2],
       [1, 1],
       [2, 0]])

pythia.index.sort_index_array(indices)[source]

Sort multiindices and remove duplicates.

Sort rows of indices by sum of multiindex and remove duplicate multiindices.

Parameters: indices (np.ndarray) – Index list before sorting.
Return type: ndarray
Returns: Sorted index array.

pythia.index.tensor_set(upper, lower=None)[source]

Create a tensor index set.

For given upper and lower bounds \(0 \leq \ell_m < u_m \in \mathbb{N}_0\) with \(m=1,\dots,M\in\mathbb{N}\), the tensor index set (n-D cube) is given by

\[\Lambda = \{ \mu\in\mathbb{N}_0^M \ \vert\ \ell_m \leq \mu_m \leq u_m \mbox{ for } m=1,\dots,M\}.\]

Parameters

upper (Sequence[int]) – Upper limit for each dimension of the tensor set. The tensor set does not include these values.
lower (Optional[Sequence[int]]) – Lower limit for each dimension of the tensor set. By default, all dimensions start with 0.

Return type

ndarray

Returns

Array with all possible multiindices in tensor set.

Examples

Create the tensor product multiindices \(\{0, 1\}\times\{0, 1\}\)

>>> pt.index.tensor_set([2, 2])
array([[0, 0],
       [0, 1],
       [1, 0],
       [1, 1]])

Create 3D univariate multiindices \(\{0\}\times\{1,\dots, 4\}\times\{0\}\)

>>> pt.index.tensor_set([1, 5, 1], [0, 1, 0])
array([[0, 1, 0],
       [0, 2, 0],
       [0, 3, 0],
       [0, 4, 0]])

Create 1D indices similar to np.arange(1, 5, dtype=int).reshape(-1, 1)

>>> pt.index.tensor_set([5], [1])
array([[1],
       [2],
       [3],
       [4]])

pythia.index.union(index_list)[source]

Build union of multiindex sets.

Given sparse index sets \(\Lambda_1, \dots, \Lambda_N\), compute \(\Lambda=\Lambda_1\cup\dots\cup\Lambda_N\).

Parameters: index_list (list of np.ndarray) – List of multiindex arrays.
Return type: ndarray
Returns: Array with all multiindices.

pythia.misc

pythia.parameter

PyThia classes containing Parameter information.

class pythia.parameter.Parameter(name, domain, distribution, mean=None, var=None, alpha=None, beta=None)[source]

Bases: object

Class used for stochasic parameters.

Parameters

name (str) – Parameter name.
domain (array_like) – Supported domain of the parameter distribution.
distribution (str) – Distribution identifier of the parameter.
mean (float, default=None) – Mean of parameter probability.
var (float, default=None) – Variance of parameter probability.
alpha (float, default=None) – Alpha value of Beta and Gamma distribution.
beta (float, default=None) – Beta value of Beta and Gamma distribution.

alpha: Union[int, float, floating, None] = None

beta: Union[int, float, floating, None] = None

distribution: str

domain: tuple[int | float | floating, int | float | floating]

mean: Union[int, float, floating, None] = None

name: str

var: Union[int, float, floating, None] = None

pythia.sampler

Sampler classes for generating in random samples and PDF evaluations.

class pythia.sampler.BetaSampler(domain, alpha, beta)[source]

Bases: Sampler

Sampler for univariate Beta distributed samples on given domain.

Parameters

domain (array_like) – Supported domain of distribution.
alpha (float) – Parameter for Beta distribution.
beta (float) – Parameter for Beta distribution.

property cov: float: (Co)Variance of the distribution.

property dimension: int: Dimension of the parameters.

grad_x_log_pdf(x)[source]

Evaluate gradient of log-PDF.

Note

Not yet implemented.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of gradient (vector valued) of log-PDF evaluated in x.

hess_x_log_pdf(x)[source]

Evaluate Hessian of log-PDF.

Note

Not yet implemented.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of Hessian (matrix valued) of log-PDF evaluated in x.

log_pdf(x)[source]

Evaluate log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of log-PDF evaluated in x.

property mass: float: Mass of the PDF.

property maximum: float

Maximum value of the PDF.

The maximum of the Beta distribution is given by

\[\begin{split}\max_{x\in[a,b]} f(x) = \begin{cases} \infty & \mbox{if } 0 < \alpha < 1 \mbox{ or } 0 < \beta < 1, \\ \frac{1}{(b-a)B(\alpha,\beta)} & \mbox{if } \alpha = 1 \mbox{ or } \beta = 1, \\ \frac{ (\alpha-1)^{\alpha-1}(\beta-1)^{\beta-1} }{ (\alpha+\beta-2)^{\alpha+\beta-2}(b-a)B(\alpha,\beta) } & \mbox{if } \alpha > 1 \mbox{ and } \beta > 1, \\ \end{cases}\end{split}\]

where \(B(\alpha,\beta)\) denotes the Beta-function.

property mean: float: Mean value of the distribution.

pdf(x)[source]

Evaluate PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of PDF evaluated in x.

sample(shape)[source]

Draw samples from distribution.

Parameters: shape (array_like) – Shape of the samples.
Return type: ndarray
Returns: Random samples of specified shape.

property std: float: Standard deviation of the distribution.

property var: float: Variance of the distribution.

class pythia.sampler.GammaSampler(domain, alpha, beta)[source]

Bases: Sampler

Sampler for univariate Gamma distributed samples on given domain.

Parameters

domain (array_like) – Supported domain of distribution.
alpha (float) – Parameter for Gamma distribution.
beta (float) – Parameter for Gamma distribution.

property cov: float: (Co)Variance of the distribution.

property dimension: int: Dimension of the parameters.

grad_x_log_pdf(x)[source]

Evaluate gradient of log-PDF.

Note

Not yet implemented.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of gradient (vector valued) of log-PDF evaluated in x.

hess_x_log_pdf(x)[source]

Evaluate Hessian of log-PDF.

Note

Not yet implemented.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of Hessian (matrix valued) of log-PDF evaluated in x.

log_pdf(x)[source]

Evaluate log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of log-PDF evaluated in x.

property mass: float: Mass of the PDF.

property maximum: float

Maximum value of the PDF.

The maximum of the Gamma distribution is given by

\[\begin{split}\max_{x\in[a,\infty)} f(x) = \begin{cases} \infty & \mbox{if } 0 < \alpha < 1 \\ \frac{\beta^\alpha}{\Gamma(\alpha)} & \mbox{if } \alpha = 1\\ \frac{\beta^\alpha}{\Gamma(\alpha)} \Bigl( \frac{\alpha-1}{\beta} \Bigr)^{\alpha-1} e^{1-\alpha} & \mbox{if } \alpha > 1\\ \end{cases}\end{split}\]

property mean: float: Mean value of the distribution.

pdf(x)[source]

Evaluate PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of PDF evaluated in x.

sample(shape)[source]

Draw samples from distribution.

Parameters: shape (array_like) – Shape of the samples.
Return type: ndarray
Returns: Random samples of specified shape.

property std: float: Standard deviation of the distribution.

property var: float: Variance of the distribution.

class pythia.sampler.NormalSampler(mean, var)[source]

Bases: Sampler

Sampler for univariate normally distributed samples.

Parameters

mean (float) – Mean of the Gaussian distribution.
var (float) – Variance of the Gaussian distribution.

property cov: float: (Co)Variance of the distribution.

property dimension: int: Dimension of the parameters.

grad_x_log_pdf(x)[source]

Evaluate gradient of log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of gradient (vector valued) of log-PDF evaluated in x.

hess_x_log_pdf(x)[source]

Evaluate Hessian of log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of Hessian (matrix valued) of log-PDF evaluated in x.

log_pdf(x)[source]

Evaluate log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of log-PDF evaluated in x.

property mass: float: Mass of the PDF.

property maximum: float: Maximum value of the PDF.

property mean: float: Mean value of the distribution.

pdf(x)[source]

Evaluate PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of PDF evaluated in x.

sample(shape)[source]

Draw samples from distribution.

Parameters: shape (array_like) – Shape of the samples.
Return type: ndarray
Returns: Random samples of specified shape.

property std: float: Standard deviation.

class pythia.sampler.ParameterSampler(params)[source]

Bases: Sampler

Product sampler of given parameters.

Parameters: params (list of pythia.parameter.Parameter) – list containing information of parameters.

property cov: ndarray: Covariance of the PDF.

property dimension: int: Dimension of the parameters.

grad_x_log_pdf(x)[source]

Evaluate gradient of log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of gradient (vector valued) of log-PDF evaluated in x.

hess_x_log_pdf(x)[source]

Evaluate Hessian of log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of Hessian (matrix valued) of log-PDF evaluated in x.

log_pdf(x)[source]

Evaluate log-PDF.

The log-PDF is given by the sum of the univariate log-PDFs.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of log-PDF evaluated in x.

property mass: float: Mass of the PDF.

property maximum: float: Maximum value of the PDF.

property mean: ndarray: Mean of the PDF.

pdf(x)[source]

Evaluate PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of PDF evaluated in x.

sample(shape)[source]

Draw samples from distribution.

Parameters: shape (array_like) – Shape of the samples.
Return type: ndarray
Returns: Random samples of specified shape.

weight(x)[source]

Weights of the parameter product PDF.

Parameters: x (np.ndarray) – Evaluation points.
Return type: ndarray
Returns: Array of uniform weights for samples.

class pythia.sampler.ProductSampler(sampler_list)[source]

Bases: Sampler

Tensor sampler for independent parameters.

Sampler for cartesian product samples of a list of (independent) univariate samplers.

Parameters: sampler_list (list of pythia.sampler.Sampler) – list of (univariate) Sampler objects.

property cov: ndarray: Covariance of the PDF.

property dimension: int: Dimension of the parameters.

grad_x_log_pdf(x)[source]

Evaluate gradient of log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of gradient (vector valued) of log-PDF evaluated in x.

hess_x_log_pdf(x)[source]

Evaluate Hessian of log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of Hessian (matrix valued) of log-PDF evaluated in x.

log_pdf(x)[source]

Evaluate log-PDF.

The log-PDF is given by the sum of the univariate log-PDFs.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of log-PDF evaluated in x.

property mass: float: Mass of the PDF.

property maximum: float: Maximum value of the PDF.

property mean: ndarray: Mean of the PDF.

pdf(x)[source]

Evaluate PDF.

The PDF is given by the product of the univariate PDFs.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of PDF evaluated in x.

sample(shape)[source]

Draw samples from distribution.

Parameters: shape (array_like) – Shape of the samples.
Return type: ndarray
Returns: Random samples of specified shape.

weight(x)[source]

Weights of the product PDF.

Parameters: x (np.ndarray) – Evaluation points.
Return type: ndarray
Returns: Array of uniform weights for samples.

class pythia.sampler.Sampler[source]

Bases: ABC

Base class for all continuous samplers.

abstract property cov: float | numpy.ndarray: (Co)Variance of the pdf.

abstract property dimension: int: Dimension of the ambient space.

domain: ndarray

abstract grad_x_log_pdf(x)[source]

Gradient of log-density of the samplers distribution.

Computes the gradient of the log-density of the samplers underlying distribution at the given points x.

Parameters: x (array_like of shape (..., D)) – list of points or single point. D is the objects dimension.
Return type: ndarray
Returns: Gradient values of the log-density at the points with shape (…, D).

abstract hess_x_log_pdf(x)[source]

Hessian of log-density of the samplers distribution.

Computes the Hessian of the log-density of the samplers underlying distribution at the given points x.

Parameters: x (array_like of shape (..., D)) – list of points or single point. D is the objects dimension.
Return type: ndarray
Returns: Hessian values of the log-density at the points with shape (…, D, D).

abstract log_pdf(x)[source]

Log-density of the samplers distribution.

Computes the log-density of the samplers underlying distribution at the given points x.

Parameters: x (array_like of shape (..., D)) – list of points or single point. D is the objects dimension.
Return type: ndarray
Returns: Log-density values at the points.

abstract property mass: float

Mass of the sampler distribution.

The integral of the sampler distribution over the domain of definition. If the density is normalised this value should be one.

abstract property maximum: float: Maximum of the pdf.

abstract property mean: float | numpy.ndarray: Mean value of the pdf.

abstract pdf(x)[source]

Density of the samplers distribution.

Computes the density of the samplers underlying distribution at the given points x.

Parameters: x (array_like of shape (..., D)) – list of points or single point. D is the objects dimension.
Return type: ndarray
Returns: Density values at the points.

abstract sample(shape)[source]

Random values in a given shape.

Create an array of the given shape and populate it with random samples from the samplers distribution.

Parameters: shape (array_like, optional) – The dimensions of the returned array, should all be positive. If no argument is given a single Python float is returned.
Return type: ndarray
Returns: Random values of specified shape.

class pythia.sampler.UniformSampler(domain)[source]

Bases: Sampler

Sampler for univariate uniformly distributed samples on given domain.

Parameters: domain (array_like) – Interval of support of distribution.

property cov: float: (Co)Variance of the distribution.

property dimension: int: Parameter dimension.

grad_x_log_pdf(x)[source]

Evaluate gradient of uniform log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of gradient (vector valued) of log-PDF evaluated in x.

hess_x_log_pdf(x)[source]

Evaluate Hessian of uniform log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of Hessian (matrix valued) of log-PDF evaluated in x.

log_pdf(x)[source]

Evaluate uniform log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of log-PDF evaluated in x.

property mass: float: Mass of the PDF.

property maximum: float: Maximum value of the PDF.

property mean: float: Mean value of the distribution.

pdf(x)[source]

Evaluate uniform PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of PDF evaluated in x.

sample(shape)[source]

Draw samples from uniform distribution.

Parameters: shape (array_like) – Shape of the samples.
Return type: ndarray
Returns: Random samples of specified shape.

property std: float: Standard deviation of the distribution.

property var: float: Variance of the distribution.

pythia.sampler.assign_sampler(param)[source]

Assign a univariate sampler to the given parameter.

Parameters: param (pythia.parameter.Parameter) –
Return type: Sampler
Returns: Univariate sampler.

pythia.sampler.constraint_sampling(sampler, constraints, shape)[source]

Draw samples according to algebraic constraints.

Draw samples from target distribution and discard samples that do not satisfy the constraints.

Parameters

sampler (Sampler) – Sampler to sample from.
constraints (list of callable) – list of functions that return True if sample point satisfies the constraint.

Return type

ndarray

Returns

Samples drawn from sampler satisfying the constraints.

Notes

The constaints may lead to a non-normalized density function.

pythia.sampler.get_maximum(f, domain, n=1000)[source]

Compute essential maximum of function by point evaluations.

Parameters

f (callable) – Function to evaluate. Needs to map from n-dim space to 1-dim space.
domain (array_like) – Domain to evaluate function on.
n (int, default=1000) – Number of function evaluations. Evaluations are done on a uniform grid in domain. Actual number of points may thus be a little greater.

Return type

float

Returns

Approximation of maximum of function f.

pythia.sampler.rejection_sampling(pdf, trial_sampler, scale, dimension, shape)[source]

Draw samples from pdf by rejection sampling.

Parameters

pdf (Callable) – Probability density the samples are generated from.
trial_sampler (Sampler) – Trial sampler proposal samples are drawn from.
scale (float) – Threshold parameter with pdf <= scale * trialSampler.pdf
dimension (int) – Dimension of the (input of the) pdf.
shape (array_like) – Shape of the samples.

Return type

ndarray

Returns

Random samples of specified shape.

pythia.least_squares_sampler

class pythia.least_squares_sampler.OWLSSampler(params, basis, tsa=True, trial_sampler=None, bulk=None)[source]

Bases: Sampler

Weighted Least-Squares sampler.

Given a stochastic variable \(y\in\Gamma\subset\mathbb{R}^M\) with \(y\sim\pi\), a set of multiindices \(\Lambda\subset\mathbb{N}_0^M\) and a finite subset \(\{P_\alpha\}_{\alpha\in\Lambda}\) of an orthonormal polynomial basis of \(L^2(\Gamma,\pi)\), the optimal weighted least-squares sampling distribution for a function \(u\in\operatorname{span}\{P_\alpha\ \vert\ \alpha\in\Lambda\}\) reads

\[\mathrm{d}\mu = w^{-1} \mathrm{d}\pi \qquad\mbox{with weight}\qquad w^{-1}(y) = \frac{1}{|\Lambda|}\sum_{\alpha\in\Lambda} |P_\alpha(y)|^2,\]

where \(\vert\Lambda\vert\) denotes the number of elements in \(\Lambda\).

Parameters

params (list of pythia.parameter.Parameter) – list of parameters.
basis (list) – list of basis functions.
tsa (bool, default=False) – Trial sampler adaptation. If True, a trial sampler is chosen on the distributions of parameters, if false a uniform trial sampler is used.
trial_sampler (pythia.sampler.Sampler, default=None) – Trial sampler for rejection sampling. If tsa is true and either trial_sampler or bulk are None, the trial sampler is chosen automatically.
bulk (float, defaul=None) – Scaling for trial sampler. If tsa is true and either trial_sampler or bulk are None, the trial sampler is chosen automatically.

Notes

To generate samples from the weighted least-squares distribution rejection sampling is used. For certain basis functions it is possible to choose a well-suited trial sampler for the rejection sampling, which can be enabled via setting tsa=True.

See also

pythia.sampler.OWLSUnivariateSampler, pythia.sampler.OWLSTensorSampler

References

The optimal weighted least-squares sampling is based on the results of Cohen & Migliorati [1]_.

property cov: ndarray: Covariance of the PDF.

property dimension: int: Dimension of the parameters.

grad_x_log_pdf(x)[source]

Evaluate gradient of log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of gradient (vector valued) of log-PDF evaluated in x.

hess_x_log_pdf(x)[source]

Evaluate Hessian of log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of Hessian (matrix valued) of log-PDF evaluated in x.

log_pdf(x)[source]

Evaluate log-PDF.

The log-PDF is given by the sum of the univariate log-PDFs.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of log-PDF evaluated in x.

property mass: float: Mass of the PDF.

property maximum: float: Maximum value of the PDF.

property mean: ndarray: Mean of the PDF.

pdf(x)[source]

Evaluate PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of PDF evaluated in x.

sample(shape)[source]

Draw samples from distribution.

Parameters: shape (array_like) – Shape of the samples.
Return type: ndarray
Returns: Random samples of specified shape.

weight(x)[source]

Weights for the PDF.

Parameters: x (array_like) – Points the weight function is evaluated in.
Return type: ndarray
Returns: weights of evaluation points x.

class pythia.least_squares_sampler.OWLSTensorSampler(params, deg, tsa=True)[source]

Bases: Sampler

Weighted least-squares sampler for tensor multivariate basis.

Given a stochastic variable \(y\in\Gamma\subset\mathbb{R}^M\) with \(y\sim\pi=\prod_{m=1}^M\pi_m\) for one dimensional densities \(\pi_m\), a tensor set of multiindices \(\Lambda=[d_1]\times\dots\times[d_M]\subset\mathbb{N}_0^M\), where \([d_m]=\{0,\dots,d_m-1\}\), and a finite subset \(\{P_\alpha\}_{\alpha\in\Lambda}\) of an orthonormal product polynomial basis of \(L^2(\Gamma,\pi)\), i.e., \(P_\alpha(y) = \prod_{m=1}^M P_{\alpha_m}(y_m)\), the optimal weighted least-squares sampling distribution for a function \(u\in\operatorname{span}\{P_\alpha\ \vert\ \alpha\in\Lambda\}\) reads

\[\mathrm{d}\mu = w^{-1} \mathrm{d}\pi \qquad\mbox{with weight}\qquad w^{-1}(y) = \prod_{m=1}^M \frac{1}{d_m} \sum_{\alpha_m\in[d_m]} |P_{\alpha_m}(y_m)|^2.\]

Parameters

params (list of pythia.parameter.Parameter) – Parameter list.
deg (list of int) – Polynomial degree of each component (same for all).
tsa (bool, default=True) – Trial sampler adaptation. If True, a trial sampler is chosen on the distributions of parameters, if false a uniform trial sampler is used.

Notes

To generate samples from the weighted least-squares distribution rejection sampling is used. For certain basis functions it is possible to choose a well-suited trial sampler for the rejection sampling, which can be enabled via setting tsa=True.

See also

pythia.sampler.OWLSUnivariateSampler

References

The optimal weighted least-squares sampling is based on the results in Cohen & Migliorati [1]_.

property cov: ndarray: Covariance of the PDF.

property dimension: int: Dimension of the parameters.

grad_x_log_pdf(x)[source]

Evaluate gradient of log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of gradient (vector valued) of log-PDF evaluated in x.

hess_x_log_pdf(x)[source]

Evaluate Hessian of log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of Hessian (matrix valued) of log-PDF evaluated in x.

log_pdf(x)[source]

Evaluate log-PDF.

The log-PDF is given by the sum of the univariate log-PDFs.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of log-PDF evaluated in x.

property mass: float: Mass of the PDF.

property maximum: float: Maximum value of the PDF.

property mean: ndarray: Mean of the PDF.

pdf(x)[source]

Evaluate PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of PDF evaluated in x.

sample(shape)[source]

Draw samples from distribution.

Parameters: shape (array_like) – Shape of the samples.
Return type: ndarray
Returns: Random samples of specified shape.

weight(x)[source]

Weights for the PDF.

Parameters: x (array_like) – Points the weight function is evaluated in.
Return type: ndarray
Returns: Weights of evaluation points x.

class pythia.least_squares_sampler.OWLSUnivariateSampler(param, deg, tsa=True)[source]

Bases: Sampler

Sampler for univariate optimally distributed samples on given domain.

Given a stochastic variable \(y\in\Gamma\subset\mathbb{R}\) with \(y\sim\pi\) and a finite subset \(\{P_j\}_{j=0}^{d-1}\) of an orthonormal polynomial basis of \(L^2(\Gamma,\pi)\), the optimal weighted least-squares sampling distribution for a function \(u\in\operatorname{span}\{P_j\ \vert\ j=0,\dots,d-1 \}\) reads

\[\mathrm{d}\mu = w^{-1} \mathrm{d}\pi \qquad\mbox{with weight}\qquad w^{-1}(y) = \frac{1}{d}\sum_{j=0}^{d-1}\vert P_j(y)\vert^2.\]

Parameters: domain (array_like) – Interval of support of distribution.

Notes

To generate samples from the weighted least-squares distribution rejection sampling is used. For certain basis functions it is possible to choose a well-suited trial sampler for the rejection sampling, which can be enabled via setting tsa=True.

See also

pythia.sampler.OWLSTensorSampler

References

The optimal weighted least-squares sampling is based on the results in Cohen & Migliorati [1]_.

1: Cohen, A. and Migliorati, G., “Optimal weighted least-squares methods”, SMAI Journal of Computational Mathematics 3, 181-203 (2017).

property cov: float: (Co)Variance of the distribution.

property dimension: int: Parameter dimension.

grad_x_log_pdf(x)[source]

Evaluate gradient of uniform log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of gradient (vector valued) of log-PDF evaluated in x.

hess_x_log_pdf(x)[source]

Evaluate Hessian of uniform log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of Hessian (matrix valued) of log-PDF evaluated in x.

log_pdf(x)[source]

Evaluate uniform log-PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of log-PDF evaluated in x.

property mass: float: Mass of the PDF.

property maximum: float: Maximum value of the PDF.

property mean: float: Mean value of the distribution.

pdf(x)[source]

Evaluate uniform PDF.

Parameters: x (array_like) – Evaluation points.
Return type: ndarray
Returns: Values of PDF evaluated in x.

sample(shape)[source]

Draw samples from weighted least-squares parameter distribution.

Parameters: shape (array_like) – Shape of the samples.
Return type: ndarray
Returns: Random samples of specified shape.

property std: float: Standard deviation of the distribution.

property var: float: Variance of the distribution.

weight(x)[source]

Weights for the pdf.

Parameters: x (np.ndarray) – Points the weight function is evaluated in.
Returns: w – Weights of evaluation points x.
Return type: array_like

class pythia.least_squares_sampler.WLSSampler(params, basis, tsa=True, trial_sampler=None, bulk=None)[source]: Bases: OWLSSampler

class pythia.least_squares_sampler.WLSTensorSampler(params, deg, tsa=True)[source]: Bases: OWLSTensorSampler

class pythia.least_squares_sampler.WLSUnivariateSampler(param, deg, tsa=True)[source]: Bases: OWLSUnivariateSampler

pythia.least_squares_sampler.owls_sampling_bound(dimension, stability_threshold=0.5, failure_probability=0.5)[source]

Minimum sample size for stable optimally weighted least squares approximation.

For a (weighted) linear least squares approximation problem in an \(m\)-dimensional linear space, let \(G\) denote the estimated Gramian matrix with respect to an \(L^2\)-orthonormal basis of this space. The approximation error is small when \(G\) is close to the idenity matrix. This function returns \(n_0 \in \mathbb{N}\) such that an optimal sample (cf. [1]_) of size \(n \ge n_0\) satisfies .. math:

\|G - I\|_2 \le \delta < 1

with a probability of at least \(1 - p\).

Parameters

dimension (int) – Dimension of the linear space (\(m \in \mathbb{N}\)).
stability_threshold (float, default=0.5) – Stability threshold that should be achieved (\(\delta \in (0, 1)\)).
failure_probability (float, default=0.5) – Failure probability with which the stability threshold is breached (\(p \in (0, 1)\)).

Returns

A sample size \(n\) ensuring that \(\|G - I\|_2\) exceeds the stability_threshold \(\delta\) with a failure_probability of at most \(p\).

Return type

int

Notes

This function is derived from formula (3.1) in [1]_ .. math:

p
= \mathbb{P}[\|I - G\|_2 > \delta]
\le 2 m \exp\left(-\frac{c_\delta n}{m}\right) ,

where \(m\) is the dimension, \(\delta\) is the stability_threshold.

References

1: A. Cohen, G. Migliorati. “Optimal weighted least-squares methods,” in The SMAI journal of computational mathematics, vol. 3, pp. 181-203, 2017