easyvvuq.sampling package¶
Submodules¶
easyvvuq.sampling.base module¶
easyvvuq.sampling.csv_sampler module¶
A CSV file based sampler.
Useful for cases where you want to evaluate a sampling plan generated by other software. Will take a CSV file with an appropriate header and will output it row by row.
easyvvuq.sampling.dataframe_sampler module¶
easyvvuq.sampling.empty module¶
easyvvuq.sampling.fd module¶
- class easyvvuq.sampling.fd.FDSampler(vary=None, distribution=None, perturbation=0.05, nominal_value=None, count=0, relative_analysis=False)[source]¶
Bases:
BaseSamplingElement
- property analysis_class¶
Return a corresponding analysis class.
- property n_samples¶
Number of samples (Ns) of PCE method. - When using pseudo-spectral projection method with tensored
quadrature: Ns = (p + 1)**d
When using pseudo-spectral projection method with sparce grid quadratue: Ns = bigO((p + 1)*log(p + 1)**(d-1))
When using regression method: Ns = 2*(p + d)!/p!*d!
Where: p is the polynomial degree and d is the number of uncertain parameters.
Ref: Eck et al. ‘A guide to uncertainty quantification and sensitivity analysis for cardiovascular applications’ [2016].
- sampler_name = 'FD_sampler'¶
easyvvuq.sampling.grid_sampler module¶
A grid sampler
Useful for e.g. hyperparameter search. The “vary” dict contains the values that must be considered per (hyper)parameter, for instance:
- vary = {“x1”: [0.0, 0.5, 0.1],
“x2 = [1, 3], “x3” = [True, False]}
The sampler will create a tensor grid using all specified 1D parameter values.
easyvvuq.sampling.mc_sampler module¶
- class easyvvuq.sampling.mc_sampler.MCSampler(vary, n_mc_samples, rule='latin_hypercube', **kwargs)[source]¶
Bases:
RandomSampler
This is a Monte Carlo sampler, used to compute the Sobol indices, mean and variance of the different QoIs.
- property analysis_class¶
Return a corresponding analysis class.
- saltelli(n_mc)[source]¶
Generates a Saltelli sampling plan of n_mc*(n_params + 2) input samples needed to compute the Sobol indices. Stored in xi_mc.
Method: A. Saltelli, Making best use of model evaluations to compute sensitivity indices, Computer Physics Communications, 2002.
- Parameters:
n_mc (the number of Monte Carlo samples per input matrix. The total)
number of samples is n_mc*(n_params + 2)
- Return type:
None.
- sampler_name = 'MC_sampler'¶
easyvvuq.sampling.mcmc module¶
- class easyvvuq.sampling.mcmc.MCMCSampler(init, q, qoi, n_chains=1, likelihood=<function MCMCSampler.<lambda>>, estimator=None)[source]¶
Bases:
BaseSamplingElement
A Metropolis-Hastings MCMC Sampler.
- Parameters:
init (dict) – Initial values for each input parameter. Of the form {‘input1’: value, …}
q (function) – A function of one argument X (dictionary) that returns the proposal distribution conditional on the X.
qoi (str) – Name of the quantity of interest
n_chains (int) – Number of MCMC chains to run in paralle.
estimator (function) – To be used with replica_col argument. Outputs an estimate of some parameter when given a sample array.
- property analysis_class¶
Returns a corresponding analysis class for this sampler.
- Return type:
class
- sampler_name = 'mcmc_sampler'¶
- update(result, invalid)[source]¶
Performs the MCMC sampling procedure on the campaign.
- Parameters:
result (pandas DataFrame) – run information from previous iteration (same as collation DataFrame)
invalid (pandas DataFrame) – invalid run information (runs that cannot be executed for some reason)
- Return type:
list of rejected run ids
easyvvuq.sampling.pce module¶
- class easyvvuq.sampling.pce.PCESampler(vary=None, distribution=None, count=0, polynomial_order=4, regression=False, rule='G', sparse=False, growth=False, relative_analysis=False, nominal_value=None)[source]¶
Bases:
BaseSamplingElement
- property analysis_class¶
Return a corresponding analysis class.
- property n_samples¶
Number of samples (Ns) of PCE method. - When using pseudo-spectral projection method with tensored
quadrature: Ns = (p + 1)**d
When using pseudo-spectral projection method with sparce grid quadratue: Ns = bigO((p + 1)*log(p + 1)**(d-1))
When using regression method: Ns = 2*(p + d)!/p!*d!
Where: p is the polynomial degree and d is the number of uncertain parameters.
Ref: Eck et al. ‘A guide to uncertainty quantification and sensitivity analysis for cardiovascular applications’ [2016].
- sampler_name = 'PCE_sampler'¶
easyvvuq.sampling.qmc module¶
This sampler is meant to be used with the QMC Analysis module.
- class easyvvuq.sampling.qmc.QMCSampler(vary, n_mc_samples, count=0)[source]¶
Bases:
BaseSamplingElement
- property analysis_class¶
Return a corresponding analysis class.
- property n_samples¶
Returns the number of samples in this sampler.
- Returns:
This computed with the formula (d + 2) * N, where d is the number
of uncertain parameters and N is the (estimated) number of samples
for the Monte Carlo method.
- sampler_name = 'QMC_sampler'¶
easyvvuq.sampling.quasirandom module¶
Summary¶
This module provides classes based on RandomSampler but modified in such a way that the output of the sampler is not random but is meant to be used in place of uniformly random number sequences. Usually this is used to cover the sampling space more “evenly” than a uniform random distribution would. Two methods are implemented:
https://en.wikipedia.org/wiki/Latin_hypercube_sampling https://en.wikipedia.org/wiki/Halton_sequence
- class easyvvuq.sampling.quasirandom.HaltonSampler(vary=None, count=0, max_num=0, analysis_class=None)[source]¶
Bases:
RandomSampler
- sampler_name = 'halton_sampler'¶
- class easyvvuq.sampling.quasirandom.LHCSampler(vary=None, count=0, max_num=0, analysis_class=None)[source]¶
Bases:
RandomSampler
- sampler_name = 'lhc_sampler'¶
easyvvuq.sampling.random module¶
- class easyvvuq.sampling.random.RandomSampler(vary=None, count=0, max_num=0, analysis_class=None)[source]¶
Bases:
BaseSamplingElement
- property analysis_class¶
Return a corresponding analysis class.
- n_samples()[source]¶
Returns the number of samples in this sampler. :rtype: if the user specifies maximum number of samples than return that, otherwise - error
- sampler_name = 'random_sampler'¶
easyvvuq.sampling.replica_sampler module¶
Replica Sampler
Summary¶
Primarily intended for sampling the same paramater values but with different random seed. Other uses may be possible. It takes a finite sampler and produces an infinite sampler from it. This infinite sampler loops through the parameters produced by the finite sampler and at each cycle adds a unique id number for that cycle to the parameter dictionary.
- class easyvvuq.sampling.replica_sampler.ReplicaSampler(sampler, replica_col='ensemble_id', seed_col=None, replicas=0)[source]¶
Bases:
BaseSamplingElement
Replica Sampler
- Parameters:
sampler (an instance of a class derived from BaseSamplingElement) – a finite sampler to loop over
replica_col (string) – a parameter name for the replica id
seed_col (string) – a parameter name for the input parameter that specifies the RNG seed
replicas (int) – number of replicas, if zero will result in an infinite sampler
- property analysis_class¶
- property inputs¶
- property iteration¶
int([x]) -> integer int(x, base=10) -> integer
Convert a number or string to an integer, or return 0 if no arguments are given. If x is a number, return x.__int__(). For floating point numbers, this truncates towards zero.
If x is not a number or if base is given, then x must be a string, bytes, or bytearray instance representing an integer literal in the given base. The literal can be preceded by ‘+’ or ‘-’ and be surrounded by whitespace. The base defaults to 10. Valid bases are 0 and 2-36. Base 0 means to interpret the base from the string as an integer literal. >>> int(‘0b100’, base=0) 4
- property qoi¶
- sampler_name = 'replica_sampler'¶
easyvvuq.sampling.sampler_of_samplers module¶
easyvvuq.sampling.simplex_stochastic_collocation module¶
THE SIMPLEX STOCASTIC COLLOCATION SAMPLER OF JEROEN WITTEVEEN (1980-2015)¶
Source:
Witteveen, J. A. S., & Iaccarino, G. (2013). Simplex stochastic collocation with ENO-type stencil selection for robust uncertainty quantification. Journal of Computational Physics, 239, 1-21.
Edeling, W. N., Dwight, R. P., & Cinnella, P. (2016). Simplex-stochastic collocation method with improved scalability. Journal of Computational Physics, 310, 301-328.
- easyvvuq.sampling.simplex_stochastic_collocation.DAFSILAS(A, b, print_message=False)[source]¶
Direct Algorithm For Solving Ill-conditioned Linear Algebraic Systems,
solves the linear system when Ax = b when A is ill conditioned.
Solves for x in the non-null subspace of the solution as described in the reference below. This method utilizes Gauss–Jordan elimination with complete pivoting to identify the null subspace of a (almost) singular matrix.
X. J. Xue, Kozaczek, K. J., Kurtzl, S. K., & Kurtz, D. S. (2000). A direct algorithm for solving ill-conditioned linear algebraic systems. Adv. X-Ray Anal, 42.
- class easyvvuq.sampling.simplex_stochastic_collocation.SSCSampler(vary=None, max_polynomial_order=4)[source]¶
Bases:
BaseSamplingElement
Simplex Stochastic Collocation sampler
- check_LEC(p_j, v, S_j, n_mc, max_jobs=4)[source]¶
Check the Local Extremum Conserving propery of all simplex elements.
- Parameters:
p_j (array, shape (n_e,)) – The polynomial order of each element.
v (array, shape (N + 1,)) – The (scalar) code outputs. #TODO:modify when vectors are allowed
S_j (array, shape (n_e, n_s)) – The indices of all nearest neighbours points of each simplex j=1,..,n_e, ordered from closest to the neighbour that furthest away. The first n_xi + 1 indeces belong to the j-th simplex itself.
n_mc (int) – The number of Monte Carlo samples to use in checking the LEC conditions.
max_jobs (int) – The number of LEC check (one per element) that can be run in parallel.
- Return type:
None.
- check_LEC_j(p_j, v, S_j, n_mc, queue)[source]¶
Check the LEC conditin of the j-th interpolation stencil.
- Parameters:
p_j (int) – The polynomial order of the j-th stencil.
v (array) – The code samples.
S_j (array, shape (N + 1,)) – The interpolation stencil of the j-th simplex element, expressed as the indices of the simplex points.
n_mc (int) – The number of Monte Carlo samples to use in checking the LEC conditions.
queue (multiprocessing queue object) – Used to store the results.
- Returns:
in the queue)
- Return type:
None, results (polynomial order and element indices are stored
- compute_ENO_stencil(p_j, S_j, el_idx, max_jobs=4)[source]¶
Compute the Essentially Non-Oscillatory stencils. The idea behind ENO stencils is to have higher degree interpolation stencils up to a thin layer of simplices containing the discontinuity. For a given simplex, its ENO stencil is created by locating all the nearest-neighbor stencils that contain element j , and subsequently selecting the one with the highest polynomial order p_j . This leads to a Delaunay triangulation which captures the discontinuity better than its nearest-neighbor counterpart.
- Parameters:
p_j (array, shape (n_e,)) – The polynomial order of each simplex element.
S_j (array, shape (n_e, n_s)) – The indices of all nearest neighbours points of each simplex j=1,..,n_e, ordered from closest to the neighbour that furthest away. The first n_xi + 1 indeces belong to the j-th simplex itself.
el_idx (dict) – The element indices for each interpolation stencil. el_idx[2] gives the elements indices of the 3rd interpolation stencil. The number of elements is determined by the local polynomial order.
max_jobs (int, optional) – The number of ENO stencils that are computed in parallel. The default is 4.
- Returns:
ENO_S_j (array, shape (n_e, n_s)) – The ENO stencils for each element.
p_j (array, shape (n_e,)) – The new polynomial order of each element.
el_idx (dict) – The new element indices for each interpolation stencil.
- compute_ENO_stencil_j(p_j, S_j, xi_centers, j, el_idx, queue)[source]¶
Compute the ENO stencil of the j-th element.
- Parameters:
p_j (array, shape (n_e,)) – The polynomial order of each simplex element.
S_j (array, shape (n_e, n_s)) – The indices of all nearest neighbours points of each simplex j=1,..,n_e, ordered from closest to the neighbour that furthest away. The first n_xi + 1 indeces belong to the j-th simplex itself.
xi_centers (array, shape (n_e, )) – The center of each simplex.
j (int) – The index of the current simplex.
el_idx (dict) – The element indices for each interpolation stencil. el_idx[2] gives the elements indices of the 3rd interpolation stencil. The number of elements is determined by the local polynomial order.
queue (multiprocessing queue object) – Used to store the results.
- Returns:
in the queue)
- Return type:
None, results (polynomial order and element indices are stored
- compute_Psi(xi_Sj, pmax)[source]¶
Compute the Vandermonde matrix Psi, given N + 1 points xi from the j-th interpolation stencil, and a multi-index set of polynomial orders |i| = i_1 + … + i_n_xi <= polynomial order.
- Parameters:
xi_Sj (array, shape (N + 1, n_xi)) – The simplex n_xi-dimensional points of the j-th interpolation stencil S_j.
pmax (int) – The max polynomial order of the local stencil.
- Returns:
Psi – The Vandermonde interpolation matrix consisting of monomials xi_1 ** i_1 + … + xi_{n_xi} ** i_{n_xi}.
- Return type:
array, shape (N + 1, N + 1)
- compute_eps_bar_j(p_j, prob_j)[source]¶
Compute the geometric refinement measure ar{eps}_j for all elements, Elements with the largest values will be selected for refinement.
- Parameters:
p_j (array, shape (n_e,)) – The polynomial order of each simplex element.
prob_j (array, shape (n_e,)) – The probability of each simplex element.
- Returns:
eps_bar_j (array, shape (n_e,)) – The geometric refinement measure.
vol_j (array, shape (n_e,)) – The volume of each simplex element.
- compute_i_norm_le_pj(p_max)[source]¶
Compute the multi-index set {i | |i| = i_1 + … + i_{n_xi} <= p}, for p = 1,…,p_max
- Parameters:
p_max (int) – The max polynomial order of the index set.
- Returns:
i_norm_le_pj – The Np + 1 multi indices per polynomial order.
- Return type:
dict
- compute_probability()[source]¶
Compute the probability Omega_j for all simplex elements;
Omega_j = int p(xi) dxi,
with integration separately over each simplex using Monte Carlo sampling.
- Returns:
prob_j – The probabilities of each simplex.
- Return type:
array, size (n_e,)
- compute_stencil_j()[source]¶
Compute the nearest neighbor stencils of all simplex elements. The distance to all points are measured with respect to the cell center of each element.
- Returns:
S_j – The indices of all nearest neighbours points of each simplex j=1,..,n_e, ordered from closest to the neighbour that furthest away. The first n_xi + 1 indeces belong to the j-th simplex itself.
- Return type:
array, shape (n_e, n_s)
- compute_sub_simplex_vertices(simplex_idx)[source]¶
Compute the vertices of the sub-simplex. The sub simplex is contained in the larger simplex. The larger simplex is refined by randomly placing a sample within the sub simplex.
- Parameters:
simplex_idx (int) – The index of the simplex
- Returns:
xi_sub_jl – The vertices of the sub simplex.
- Return type:
array, size (n_xi + 1, n_xi)
- compute_surplus_k(xi_k_jref, S_j, p_j, v, v_k_jref)[source]¶
Compute the hierachical surplus at xi_k_jref (the refinement location), defined as the difference between the new code sample and the (old) surrogate prediction at the refinement location.
- Parameters:
xi_k_jref (array, shape (n_xi,)) – The refinement location.
S_j (array, shape (n_e, n_s)) – The indices of all nearest neighbours points of each simplex j=1,..,n_e, ordered from closest to the neighbour that furthest away. The first n_xi + 1 indeces belong to the j-th simplex itself.
p_j (array, shape (n_e,)) – The polynomial order of each simplex element.
v (array, shape (N + 1,)) – The (scalar) code outputs. #TODO:modify when vectors are allowed
v_k_jref (float #TODO:modify when vectors are allowed) – The code prediction at the refinement location.
- Returns:
surplus – The hierarchical suplus
- Return type:
float #TODO:modify when vectors are allowed
- compute_vol()[source]¶
Use analytic formula to compute the volume of all simplices.
https://en.wikipedia.org/wiki/Simplex#Volume
- Returns:
vols – The volumes of the n_e simplices.
- Return type:
array, size (n_e,)
- compute_xi_center_j()[source]¶
Compute the center of all simplex elements.
- Returns:
xi_center_j – The cell centers of the n_e simplices.
- Return type:
array, size (n_e,)
- find_boundary_simplices()[source]¶
Find the simplices that are on the boundary of the hypercube .
- Returns:
idx – Indices of the boundary simplices.
- Return type:
array
- find_pmax(n_s)[source]¶
Finds the maximum polynomial stencil order that can be sustained given the current number of code evaluations. The stencil size required for polynonial order p is given by;
stencil size = (n_xi + p)! / (n_xi!p!)
where n_xi is the number of random inputs. This formula is used to find p.
- Parameters:
n_xi (int) – Number of random parameters.
n_s (int) – Number of code samples.
- Returns:
The highest polynomial order that can be sustained given the number of code evaluations.
- Return type:
int
- find_simplices(S_j)[source]¶
Find the simplex element indices that are in point stencil S_j.
- Parameters:
S_j (array, shape (N + 1,)) – The interpolation stencil of the j-th simplex element, expressed as the indices of the simplex points.
- Returns:
idx – The element indices of stencil S_j.
- Return type:
array
- init_grid()[source]¶
Create the initial n_xi-dimensional Delaunay discretization
- Returns:
tri – Initial triagulation of 2 ** n_xi + 1 points.
- Return type:
scipy.spatial.qhull.Delaunay
- load_state(filename)[source]¶
Load the state of the sampler from a pickle file.
- Parameters:
filename (string) – File name.
- Return type:
None.
- property n_dimensions¶
Returns the number of uncertain random variables
- property n_elements¶
Returns the number of simplex elements
- property n_samples¶
Returns the number of samples code samples.
- sample_inputs(n_mc)[source]¶
Draw n_mc random values from the input distribution.
- Parameters:
n_mc (int) – The number of Monte Carlo samples.
- Returns:
xi – n_mc random samples from the n_xi input distributions.
- Return type:
array, shape(n_mc, n_xi)
- sample_simplex(n_mc, xi_k_jl, check=False)[source]¶
Use an analytical map from n_mc uniformly distributed points in the n_xi-dimensional hypercube, to uniformly distributed points in the target simplex described by the nodes xi_k_jl.
Derivation: Edeling, W. N., Dwight, R. P., & Cinnella, P. (2016). Simplex-stochastic collocation method with improved scalability. Journal of Computational Physics, 310, 301-328.
- Parameters:
n_mc (int) – The number of Monte Carlo samples.
xi_k_jl (array, shape (n_xi + 1, n_xi)) – The nodes of the target simplex.
check (bool, optional) – Check is the random samples actually all lie insiide the target simplex. The default is False.
- Returns:
P – Uniformly distributed points inside the target simplex.
- Return type:
array, shape (n_mc, n_xi)
- sample_simplex_edge(simplex_idx, refined_edges)[source]¶
Refine the longest edge of a simplex.
# TODO: is allright for 2D, but does this make sense in higher dims?
- Parameters:
simplex_idx (int) – The index of the simplex.
refined_edges (list) – Contains the pairs of the point indices, corresponding to edges that have been refined in the current iteration. Simplices share edges, and this list is used to prevent refining the same edge twice within the same iteration of the SSC algorihm.
- Returns:
xi_new (array, shape (n_xi,)) – The newly added point (along the longest edge).
refined_edges (list) – The updated refined_edges list.
already_refined (bool) – Returns True if the edge already has been refined.
- sampler_name = 'ssc_sampler'¶
- save_state(filename)[source]¶
Save the state of the sampler to a pickle file.
- Parameters:
filename (string) – File name.
- Return type:
None.
- set_pmax_cutoff(pmax_cutoff)[source]¶
Set the maximum allowed polynomial value of the surrogate.
- Parameters:
p_max_cutoff (int) – The max polynomial order.
- Return type:
None.
- surrogate(xi, S_j, p_j, v)[source]¶
Evaluate the SSC surrogate at xi.
- Parameters:
xi (array, shape (n_xi,)) – The location in the input space at which to evaluate the surrogate.
S_j (array, shape (n_e, n_s)) – The indices of all nearest neighbours points of each simplex j=1,..,n_e, ordered from closest to the neighbour that furthest away. The first n_xi + 1 indeces belong to the j-th simplex itself.
p_j (array, shape (n_e,)) – The polynomial order of each simplex element.
v (array, shape (N + 1,)) – The (scalar) code outputs. #TODO:modify when vectors are allowed
- Return type:
None.
- update_Delaunay(new_points)[source]¶
Update the Delaunay triangulation with P new points.
- Parameters:
new_points (array, shape (P, n_xi)) – P new n_xi-dimensional points.
- Return type:
None.
- w_j(xi, c_jl, pmax)[source]¶
Compute the surrogate local interpolation at point xi.
# TODO: right now this assumes a scalar output. Modify the code for vector-valued outputs.
- Parameters:
xi (array, shape (n_xi,)) – A point inside the stochastic input domain.
c_jl (array, shape (N + 1,)) – The interpolation coefficients of the j-th stencil, with l = 0, …, N.
pmax (int) – The max polynomial order of the local stencil.
- Returns:
w_j_at_xi – The surrogate prediction at xi.
- Return type:
float
easyvvuq.sampling.stochastic_collocation module¶
- class easyvvuq.sampling.stochastic_collocation.SCSampler(vary=None, polynomial_order=4, quadrature_rule='G', count=0, growth=False, sparse=False, midpoint_level1=True, dimension_adaptive=False)[source]¶
Bases:
BaseSamplingElement
Stochastic Collocation sampler
- property analysis_class¶
Return a corresponding analysis class.
- check_max_quad_level()[source]¶
If a discrete variable is specified, there is the possibility of non unique collocation points if the quadrature order is high enough. This subroutine prevents that.
NOTE: Only detects cp.DiscreteUniform thus far
The max quad orders are stores in self.max_quad_order
- Return type:
None
- compute_1D_points_weights(L, N)[source]¶
Computes 1D collocation points and quad weights, and stores this in self.xi_1d, self.wi_1d.
- Parameters:
L ((int) the max polynomial order of the (sparse) grid)
N ((int) the number of uncertain parameters)
- Return type:
None.
- compute_sparse_multi_idx(L, N)[source]¶
computes all N dimensional multi-indices l = (l1,…,lN) such that |l| <= L + N - 1, i.e. a simplex set: 3 * 2 * * (L=3 and N=2) 1 * * *
1 2 3
Here |l| is the internal sum of i (l1+…+lN)
- look_ahead(current_multi_idx)[source]¶
The look-ahead step in (dimension-adaptive) sparse grid sampling. Allows for anisotropic sampling plans.
Computes the admissible forward neighbors with respect to the current level multi-indices. The admissible level indices l are added to self.admissible_idx. The code will be evaluated next iteration at the new collocation points corresponding to the levels in admissble_idx.
Source: Gerstner, Griebel, “Numerical integration using sparse grids”
- Parameters:
current_multi_idx (array of the levels in the current iteration)
of the sparse grid.
- Return type:
None.
- property n_samples¶
Number of samples (Ns) of SC method. - When using tensor quadrature: Ns = (p + 1)**d - When using sparid: Ns = bigO((p + 1)*log(p + 1)**(d-1)) Where: p is the polynomial degree and d is the number of uncertain parameters.
Ref: Eck et al. ‘A guide to uncertainty quantification and sensitivity analysis for cardiovascular applications’ [2016].
- next_level_sparse_grid()[source]¶
Adds the points of the next level for isotropic hierarchical sparse grids.
- Return type:
None.
- sampler_name = 'sc_sampler'¶
easyvvuq.sampling.sweep module¶
easyvvuq.sampling.transformations module¶
Module contents¶
Classes implementing the sampling element for EasyVVUQ
Summary¶
Samplers in the context of EasyVVUQ are classes that generate sequences of parameter dictionaries. These dictionaries are then used to create input files for the simulations.