Python API
Calibration toolbox
This module contains the Bayesian calibration class.
- class BayesianCalibration(system, calibration, num_iter, curr_iter, save_fig)
Bases:
object
This is the Bayesian calibration class.
A Bayesian calibration class consists of the inference class and the (re)sampling class. For instance, in GrainLearning, we have a calibration method called “iterative Bayesian filter” which consists of “sequential Monte Carlo” for model parameter estimation and “variational Gaussian mixture” for resampling.
There are two ways of initializing a calibration toolbox class.
Method 1 - dictionary style (recommended)
bayesian_calibration = BayesianCalibration.from_dict( { "num_iter": 8, "system": { "system_type": DynamicSystem, "model_name": "test", "param_names": ["a", "b"], "param_min": [0, 0], "param_max": [1, 10], "num_samples": 10, "obs_data": [2,4,8,16], "ctrl_data": [1,2,3,4], "callback": run_sim, }, "calibration": { "inference": {"ess_target": 0.3}, "sampling": {"max_num_components": 1}, }, "save_fig": -1, } )
or
Method 2 - class style
bayesian_calibration = BayesianCalibration( num_iter = 10, system = DynamicSystem(...), calibration = IterativeBayesianFilter(...) save_fig = -1 )
- Parameters:
system (
Type
[DynamicSystem
]) – A dynamic system whose observables and hidden states evolve dynamically over “time”calibration (
Type
[IterativeBayesianFilter
]) – An iterative Bayesian Filter that iteratively sample the parameter spacenum_iter (
int
) – Number of iteration stepscurr_iter (
int
) – Current iteration stepsave_fig (
int
) – Flag for skipping (-1), showing (0), or saving (1) the figures
- classmethod from_dict(obj)
An alternative constructor to allow choosing a system type (e.g., dynamic system or IO dynamic system)
- Parameters:
obj (
Dict
) – a dictionary containing the keys and values to construct a BayesianCalibration object- Returns:
a BayesianCalibration object
- get_most_prob_params()
Return the most probable set of parameters
- Returns:
Estimated parameter values
- load_and_process(sigma=0.1)
Load existing simulation data and compute the posterior distribution using an assumed sigma
- Parameters:
sigma (
float
) – assumed uncertainty coefficient, defaults to 0.1
- load_and_run_one_iteration()
Load existing simulation data and run Bayesian calibration for one iteration Note the maximum uncertainty sigma_max is solved to reach a certain effective sample size ess_target, unlike being assumed as an input for load_and_process(…)
- load_system()
Load existing simulation data from disk into the dynamic system
- plot_uq_in_time()
Plot the evolution of uncertainty moments and distribution over time
- resample()
Learn and resample from a proposal distribution todo this should be refactored
- Returns:
Combinations of resampled parameter values
- run()
This is the main calibration loop which does the following steps 1. First iteration of Bayesian calibration starts with a Halton sequence 2. Iterations continue by resampling the parameter space until certain criteria are met.
- run_one_iteration(index=-1)
Run Bayesian calibration for one iteration.
- Parameters:
index (
int
) – iteration step, defaults to -1
-
calibration:
Type
[IterativeBayesianFilter
] Calibration method (e.g, Iterative Bayesian Filter)
-
curr_iter:
int
= 0 Current calibration step
-
num_iter:
int
Number of iterations
-
save_fig:
int
= -1 Flag to save figures
-
system:
Type
[DynamicSystem
] Dynamic system whose parameters or hidden states are being inferred
Dynamic systems
This module contains various classes of state-space systems that describe the time evolution of the system’s hidden state based on partial observations from the real world.
- class DynamicSystem(obs_data, num_samples, param_min, param_max, ctrl_data=None, obs_names=None, ctrl_name=None, inv_obs_weight=None, sim_name=None, sim_data=None, callback=None, param_data=None, param_names=None, sigma_max=1000000.0, sigma_tol=0.001)
Bases:
object
This is the dynamical system class.
A dynamical system (also known as a state-space system) describes the time evolution of the system’s (hidden) state and observation using the following equations:
\[ \begin{align}\begin{aligned}x_t & = f(x_{t−1}) + q_{t−1}\\y_t & = h(x_t) + r_t\end{aligned}\end{align} \]where \(x_t\) is the hidden state, \(y_t\) is the observation, both represented as random processes, \(f\) is the state transition function, \(h\) is the observation function, \(q_{t−1}\) is the process noise, and \(r_t\) is the observation noise.
In the context of Bayesian parameter estimation, \(f\) is the model that describe the physical process and \(h\) is the model that describe the relationship between the hidden state and observation. In the simplest case, \(h\) is an identity matrix which indicate the one-to-one relationship between \(x_t\) and \(y_t\).
Therefore, the
DynamicSystem
class is used to encapsulate the observation data and the simulation data, which require specifying the number of samples, the lower and upper bound of the parameters, and the callback function that runs the forward predictions.There are two ways of initializing the class.
Method 1 - dictionary style
system_cls = DynamicSystem.from_dict( { "param_names": ["a", "b"], "param_min": [0, 0], "param_max": [1, 10], "num_samples": 14, "obs_data": y_obs, "ctrl_data": y_ctrl, "callback": run_sim } )
or
Method 2 - class style
system_cls = DynamicSystem( param_names = ["a", "b"], param_min = [0, 0], param_max = [1, 10], num_samples = 14, obs_data = y_obs, ctrl_data = y_ctrl, callback = run_sim )
You can pass the simulation data to the class using the
DynamicSystem.set_sim_data()
method.system_cls.set_sim_data(x)
The simulation data is a numpy array of shape (num_samples, num_obs, num_steps).
- Parameters:
obs_data (
ndarray
) – observation or reference datanum_samples (
int
) – Sample size of the ensemble of model evaluationsparam_min (
List
[float
]) – List of parameter lower boundsparam_max (
List
[float
]) – List of parameter Upper boundscallback (
Optional
[Callable
]) – Callback function, defaults to Nonectrl_data (
Optional
[ndarray
]) – control data (e.g, time), defaults to None, optionalobs_names (
Optional
[List
[str
]]) – Column names of the observation data, defaults to None, optionalctrl_name (
Optional
[str
]) – Column name of the control data, defaults to None, optionalinv_obs_weight (
Optional
[List
[float
]]) – Inverse of the observation weight, defaults to None, optionalparam_data (
Optional
[ndarray
]) – Parameter data, defaults to None, optionalparam_names (
Optional
[List
[str
]]) – Parameter names, defaults to None, optionalsim_data (
Optional
[ndarray
]) – Simulation data, defaults to None, optionalsigma_max (
float
) – Maximum uncertainty, defaults to 1.0e6, optionalsigma_tol (
float
) – Tolerance of the estimated uncertainty, defaults to 1.0e-3, optionalsim_name (
Optional
[str
]) – Name of the simulation, defaults to ‘sim’, optional
- classmethod backup_sim_data()
Virtual function to backup simulation data
- compute_estimated_params(posteriors)
Compute the estimated means and uncertainties for the parameters.
This function is vectorized for all time steps
- Parameters:
posteriors (
array
) – Posterior distribution of shape (num_steps, num_samples)
- compute_inv_normalized_sigma()
Compute the inverse of the matrix that apply different weights on the observables
- classmethod from_dict(obj)
Initialize the class using a dictionary style
- Parameters:
obj (
dict
) – Dictionary object- Returns:
DynamicSystem: DynamicSystem object
- get_inv_normalized_sigma()
Get the inverse of the matrix that apply different weights on the observables
- classmethod get_sim_data_files(curr_iter)
Virtual function to get simulation data files from disk
- Parameters:
curr_iter (int) –
- classmethod load_param_data(curr_iter)
Virtual function to load param data from disk
- Parameters:
curr_iter (int) –
- classmethod load_sim_data()
Virtual function to load simulation data
- reset_inv_normalized_sigma()
Reset the inverse of the weighting matrix to None
- run(**kwargs)
Run the callback function
TODO design a better wrapper to avoid kwargs? :type kwargs: :param kwargs: keyword arguments to pass to the callback function
- set_sim_data(data)
Set the simulation data
- Parameters:
data (
list
) – simulation data of shape (num_samples, num_obs, num_steps)
- classmethod write_params_to_table(curr_iter)
Write the parameter data into a text file
-
callback:
Callable
Callback function to run the forward predictions
-
ctrl_data:
ndarray
Control dataset (num_ctrl, num_steps)
-
ctrl_name:
str
Observation control (e.g., time)
-
estimated_params:
ndarray
= None Estimated parameter as the first moment of the distribution (\(x_\mu = \sum_i w_i * x_i\))
-
estimated_params_cv:
ndarray
= None Parameter uncertainty as the second moment of the distribution (\(x_\sigma = \sum_i w_i \cdot (x_i - x_\mu)^2 / x_\mu\))
-
inv_obs_weight:
List
[float
] Inverse of the observation weight
-
num_ctrl:
int
Number of control data (identical between simulation and observation)
-
num_obs:
int
Number of observations in the dataset
-
num_params:
int
Number of unknown parameters
-
num_samples:
int
Number of samples (usually specified by user)
-
num_steps:
int
Number of steps or sequence size in the dataset
-
obs_data:
ndarray
Observation (or reference) data of shape (num_obs, num_steps)
-
obs_names:
List
[str
] Observation keys
-
param_data:
ndarray
Parameter data of shape (num_samples, num_params)
-
param_data_prev:
ndarray
Parameter data of previous iteration
-
param_max:
List
Upper bound of the parameters
-
param_min:
List
Lower bound of the parameters
-
param_names:
List
[str
] Names of the parameters
-
sigma_max:
float
= 1000000.0 Maximum value of the uncertainty
-
sigma_min:
float
= 1e-06 Minimum value of the uncertainty
-
sigma_tol:
float
= 0.001 Sigma tolerance
-
sim_data:
ndarray
Simulation data of shape (num_samples, num_obs, num_steps)
-
sim_name:
str
= 'sim' Name of the simulation (e.g., sim)
- class IODynamicSystem(sim_name, sim_data_dir, sim_data_file_ext, obs_data_file, obs_names, ctrl_name, num_samples, param_min, param_max, obs_data=None, ctrl_data=None, inv_obs_weight=None, sim_data=None, callback=None, param_data_file=None, param_data=None, param_names=None)
Bases:
DynamicSystem
This is the I/O dynamic system class derived from the dynamic system class. Extra functionalities are added to handle I/O operations.
There are two ways of initializing the class.
Method 1 - dictionary style
system_cls = IODynamicSystem.from_dict( { "system_type": IODynamicSystem, "param_min": [0, 0], "param_max": [1, 10], "param_names": ['a', 'b'], "num_samples": 14, "obs_data_file": 'obs_data.txt', "obs_names": ['y_obs'], "ctrl_name": 'y_ctrl "sim_name": 'linear', "sim_data_file_ext": '.txt', "callback": run_sim } )
or
Method 2 - class style
system_cls = IODynamicSystem( param_min = [0, 0], param_max = [1, 10], param_names = ['a', 'b'], num_samples = 14, obs_data_file = 'obs_data.txt', obs_names = ['y_obs'], ctrl_name = 'y_ctrl', sim_name = 'linear', sim_data_file_ext = '.txt', callback = run_sim )
- Parameters:
param_min (
List
[float
]) – List of parameter lower boundsparam_max (
List
[float
]) – List of parameter Upper boundsparam_names (
Optional
[List
[str
]]) – Parameter names, defaults to Nonenum_samples (
int
) – Sample size of the ensemble of model evaluationsobs_data_file (
str
) – Observation data file, defaults to Noneobs_names (
List
[str
]) – Column names of the observation data, defaults to Nonectrl_name (
str
) – Column name of the control data, defaults to Nonesim_name (
str
) – Name of the simulation, defaults to ‘sim’sim_data_dir (
str
) – Simulation data directory, defaults to ‘./sim_data’sim_data_file_ext (
str
) – Simulation data file extension, defaults to ‘.npy’callback (
Optional
[Callable
]) – Callback function, defaults to Noneparam_data_file (
Optional
[str
]) – Parameter data file, defaults to None, optionalobs_data (
Optional
[ndarray
]) – observation or reference data, optionalctrl_data (
Optional
[ndarray
]) – Control data (e.g, time), defaults to None, optionalinv_obs_weight (
Optional
[List
[float
]]) – Inverse of the observation weight, defaults to None, optionalparam_data (
Optional
[ndarray
]) – Parameter data, defaults to None, optionalsim_data (
Optional
[ndarray
]) – Simulation data, defaults to None, optional
- backup_sim_data()
Backup simulation data files to a backup directory.
- compute_estimated_params(posteriors)
Compute the estimated means and uncertainties for the parameters.
This function is vectorized for all time steps
- Parameters:
posteriors (
array
) – Posterior distribution of shape (num_steps, num_samples)
- compute_inv_normalized_sigma()
Compute the inverse of the matrix that apply different weights on the observables
- classmethod from_dict(obj)
Initialize the class using a dictionary style
- Parameters:
obj (
dict
) – Dictionary object- Return IODynamicSystem:
IODynamicSystem object
- get_inv_normalized_sigma()
Get the inverse of the matrix that apply different weights on the observables
- get_obs_data()
Get the observation data from the observation data file.
Separate the control data from the observation data if the name of control variable is given. Otherwise, the observation data is the entire data in the observation data file.
- get_sim_data_files(curr_iter=0)
Get the simulation data files from the simulation data directory.
- Parameters:
curr_iter (
int
) – Current iteration number, default to 0.
- load_param_data(curr_iter=0)
Load parameter data from a table written in a text file.
- Parameters:
curr_iter (
int
) – Current iteration number, default to 0.
- load_sim_data()
Load the simulation data from the simulation data files.
The function does the following: 1. Load simulation data into an IO dynamic system object 2. Check if parameter values read from the table matches those used to creat the simulation data
- reset_inv_normalized_sigma()
Reset the inverse of the weighting matrix to None
- run(**kwargs)
Run the callback function
TODO design a better wrapper to avoid kwargs? :type kwargs: :param kwargs: keyword arguments to pass to the callback function
- set_sim_data(data)
Set the simulation data
- Parameters:
data (
list
) – simulation data of shape (num_samples, num_obs, num_steps)
- write_params_to_table(curr_iter)
Write the parameter data into a text file.
- Parameters:
curr_iter (
int
) – Current iteration number, default to 0.- Return param_data_file:
The name of the parameter data file
-
callback:
Callable
Callback function to run the forward predictions
-
ctrl_data:
ndarray
Control dataset (num_ctrl, num_steps)
-
ctrl_name:
str
Observation control (e.g., time)
-
estimated_params:
ndarray
= None Estimated parameter as the first moment of the distribution (\(x_\mu = \sum_i w_i * x_i\))
-
estimated_params_cv:
ndarray
= None Parameter uncertainty as the second moment of the distribution (\(x_\sigma = \sum_i w_i \cdot (x_i - x_\mu)^2 / x_\mu\))
-
inv_obs_weight:
List
[float
] Inverse of the observation weight
-
num_ctrl:
int
Number of control data (identical between simulation and observation)
-
num_obs:
int
Number of observations in the dataset
-
num_params:
int
Number of unknown parameters
-
num_samples:
int
Number of samples (usually specified by user)
-
num_steps:
int
Number of steps or sequence size in the dataset
-
obs_data:
ndarray
Observation (or reference) data of shape (num_obs, num_steps)
-
obs_data_file:
str
Name of the observation data file
-
obs_names:
List
[str
] Observation keys
-
param_data:
ndarray
Parameter data of shape (num_samples, num_params)
-
param_data_file:
str
Name of the parameter data file
-
param_data_prev:
ndarray
Parameter data of previous iteration
-
param_max:
List
Upper bound of the parameters
-
param_min:
List
Lower bound of the parameters
-
param_names:
List
[str
] Names of the parameters
-
sigma_max:
float
= 1000000.0 Maximum value of the uncertainty
-
sigma_min:
float
= 1e-06 Minimum value of the uncertainty
-
sigma_tol:
float
= 0.001 Sigma tolerance
-
sim_data:
ndarray
Simulation data of shape (num_samples, num_obs, num_steps)
-
sim_data_dir:
str
= './sim_data' Name of the directory where simulation data is stored
-
sim_data_file_ext:
str
= '.npy' Extension of simulation data files
-
sim_data_files:
List
[str
] Simulation data files (num_samples)
-
sim_name:
str
= 'sim' Name of the simulation (e.g., sim)
Iterative Bayesian Filter
This module contains various Bayesian filtering classes by mixing the inference and sampling methods
from the inference
and sampling
modules.
- class IterativeBayesianFilter(inference, sampling, ess_tol=0.01, initial_sampling='halton', proposal=None, proposal_data_file=None)
Bases:
object
This is the Iterative Bayesian Filter class.
The idea is to solve an inverse problem repeatedly by refining the previous knowledge of the dynamical system. An iterative Bayesian filter is essentially a sequential Monte Carlo filter that iteratively samples the parameter space, using a proposal distribution trained with non-Gaussian mixture models. It consists of the inference class and the (re)sampling class. The inference class is a sequential Monte Carlo filter, and the resampling class is a Gaussian mixture model. The Gaussian mixture model is trained with the previous ensemble (i.e., samples and associated weights).
The steps of the iterative Bayesian filter are as follows:
Initialize
the ensemble of model (including parameter) states using a low-discrepancy sequence.Compute the posterior distribution
of model states such that the target effective sample size is reached.Generate new samples
from a proposal density that is trained with the previous ensemble (i.e., samples and associated weights).Run the model with the new samples and update the ensemble of model states using
DynamicSystem.run()
andBayesianCalibration.load_system()
.Repeat steps 2-4 until convergence. The iterations are done in the Bayesian Calibration class by calling
BayesianCalibration.run_one_iteration()
.
There are two ways of initializing the class.
Method 1 - dictionary style
ibf_cls = IterativeBayesianFilter.from_dict( { "inference":{ "ess_target": 0.3, "scale_cov_with_max": True }, "sampling":{ "max_num_components": 2 } } )
or
Method 2 - class style
system_cls = IterativeBayesianFilter( inference = SMC(...), sampling = GaussianMixtureModel(...) )
- Parameters:
inference (
Type
[SMC
]) – Sequential Monte Carlo class (SMC)sampling (
Type
[GaussianMixtureModel
]) – Gaussian Mixture Model class (GMM)initial_sampling (
str
) – The initial sampling method, defaults to Haltoness_tol (
float
) – Tolerance for the target effective sample size to converge, defaults to 1.0e-2proposal (
Optional
[ndarray
]) – A proposal distribution to sample the state-parameter space, defaults to Noneproposal_data_file (
Optional
[str
]) – Pickle that stores the previously trained proposal distribution, defaults to None
- add_curr_param_data_to_list(param_data)
Add the current parameter samples to the list of parameter samples.
- Parameters:
param_data (
ndarray
) – Current parameter samples of shape (num_samples, num_params)
- classmethod from_dict(obj)
Initialize the class using a dictionary style
- Parameters:
obj (
dict
) – a dictionary containing the keys and values to construct an iBF object- Returns:
an iBF object
- initialize(system)
Initialize the ensemble of model (including parameter) states using a low-discrepancy sequence.
- Parameters:
system (
Type
[DynamicSystem
]) – Dynamic system class
- load_proposal_from_file(system, curr_iter)
Load the proposal density from a file.
- Parameters:
system (
Type
[IODynamicSystem
]) – Dynamic system classcurr_iter (
int
) – Current iteration number
- run_inference(system, curr_iter=0)
Compute the posterior distribution of model states such that the target effective sample size is reached.
- Parameters:
system (
Type
[DynamicSystem
]) – Dynamic system classcurr_iter (
int
) – Current iteration number
- run_sampling(system)
Generate new samples from a proposal density.
The proposal density can be trained with the previous ensemble (i.e., samples and associated weights).
- Parameters:
system (
Type
[DynamicSystem
]) – Dynamic system class
- save_proposal_to_file(system, curr_iter)
Save the proposal density to a file.
- Parameters:
system (
Type
[IODynamicSystem
]) – Dynamic system classcurr_iter (
int
) – Current iteration number
- solve(system)
Run both inference and sampling for a dynamic system
The iterations are done at the level of the Bayesian Calibration class by calling
BayesianCalibration.run_one_iteration()
.- Parameters:
system (
Type
[DynamicSystem
]) – Dynamic system class
-
ess_tol:
float
= 0.01 This a tolerance to which the optimization algorithm converges.
- inference
The inference class quantify the evolution of the posterior distribution of model states over time
alias of
Type
[SMC
]
-
initial_sampling:
str
= 'halton' The non-informative distribution to draw the initial samples
-
param_data_list:
List
= [] List of the parameter samples of shape (num_samples, num_params) generated in all iterations
-
posterior:
ndarray
= None The next proposal distribution
-
proposal:
ndarray
= None The current proposal distribution
-
proposal_data_file:
str
= None The name of the file that stores the current proposal distribution
- sampling
The sampling class generates new samples from the proposal density
alias of
Type
[GaussianMixtureModel
]
-
sigma_list:
List
= [] List of sigma values optimized to satisfy the target effective sample size in all iterations
Inference
This module contains various methods for performing statistical and Bayesian inference
- class SMC(ess_target, scale_cov_with_max=True, cov_matrices=None)
Bases:
object
This is the Sequential Monte Carlo class that recursively update model states based on sequential observations using the Bayes’ theorem
There are two ways of initializing the class.
Method 1 - dictionary style
system_cls = SMC.from_dict( { "ess_target": 0.3, "scale_cov_with_max": True } )
or
Method 2 - class style
system_cls = SMC( ess_target = 0.3, scale_cov_with_max = True )
- Parameters:
ess_target (
float
) – Target effective sample size (w_0 / sum_i w_i^2) where w_0 = 1 / N_p, defaults to 0.3scale_cov_with_max (
bool
) – True if the covariance matrix is scaled with the maxima of the observations, defaults to Truecov_matrices (
Optional
[array
]) – Covariance matrices of shape (num_steps, num_obs, num_obs), defaults to None, Optional
- compute_effective_sample_size()
Compute the effective sample size
- data_assimilation_loop(sigma, system, proposal=None)
Perform data assimilation loop
- Parameters:
sigma (
float
) – Uncertaintyproposal (
Optional
[ndarray
]) – Proposal distribution from which the samples are sampled from, defaults to None, optionalsystem (
Type
[DynamicSystem
]) – Dynamic system class
- Returns:
Result of the objective function which converges to a user defined effective sample size
- classmethod from_dict(obj)
Initialize the class using a dictionary style
- Parameters:
obj (
dict
) – a dictionary containing the keys and values to construct an SMC object- Returns:
an SMC object
- get_covariance_matrices(sigma, system)
Compute the (diagonal) covariance matrices from an uncertainty that is either assumed by the user or updated by SMC to satisfy the target effective sample size.
This function is vectorized for all time steps
- Parameters:
sigma (
float
) – Uncertaintysystem (
Type
[DynamicSystem
]) – Dynamic system
- Returns:
Covariance matrices for all time steps
- static get_likelihoods(system, cov_matrices)
Compute the likelihood distributions of simulation data as a multivariate normal centered around the observation.
This function is vectorized for all time steps
- Parameters:
system (
Type
[DynamicSystem
]) – Dynamic system classcov_matrices (
array
) – Covariance matrices of shape (num_steps, num_obs, num_obs)
- Returns:
Likelihood matrices of shape (num_steps, num_samples) considering all time steps
- get_posterior_at_time(time_step=-1)
Get the posterior distribution at a certain time step
- Parameters:
time_step (
int
) – input time step, defaults to -1 (the last step in time), optional- Returns:
Posterior distribution at a certain time step
- static get_posteriors(system, likelihoods, proposal=None)
Compute the posterior distributions from the likelihood for all the time steps
This function is vectorized for all time steps
- Parameters:
system (
Type
[DynamicSystem
]) – Dynamic system classlikelihoods (
array
) – Likelihood matrices of shape (num_steps, num_samples)proposal (
Optional
[array
]) – Proposal distribution at the initial time step, defaults to None, optional
- Returns:
Posterior matrices of shape (num_steps, num_samples) considering all time steps
-
cov_matrices:
array
Covariance matrices of shape (num_steps, num_obs, num_obs)
-
ess:
array
Time evolution of the effective sample size
-
ess_target:
float
Target effective sample size
-
likelihoods:
array
Likelihood distributions of shape (num_steps, num_samples)
-
posteriors:
array
Posterior distributions of shape (num_steps, num_samples)
-
scale_cov_with_max:
bool
= True True if the covariance matrix is scaled with the maximum of the observations, defaults to True
Sampling
This module contains various methods to sample the state-parameter space of a dynamic system.
- class GaussianMixtureModel(max_num_components, weight_concentration_prior=None, covariance_type='tied', n_init=1, tol=1e-05, max_iter=100, random_state=None, init_params='kmeans', warm_start=False, expand_factor=10, slice_sampling=False)
Bases:
object
This is the wrapper class of the Bayesian Gaussian Mixture. from scikit-learn.
The BayesianGaussianMixture. class is extended to include the following functionalities:
GaussianMixtureModel.expand_and_normalize_weighted_samples()
: Converting an ensemble of samples with importance weights into an ensemble of samples with equal weightsGaussianMixtureModel.regenerate_params()
: Regenerating new samples by training and sampling from a Gaussian mixture model.
There are two ways of initializing the class.
Method 1 - dictionary style
system_cls = GaussianMixtureModel.from_dict( { "max_num_components": 2, "covariance_type": "full", "random_state": 0, slice_sampling: False, } )
or
Method 2 - class style
system_cls = GaussianMixtureModel( max_num_components = 2, covariance_type = "full", random_state = 0, slice_sampling = False, )
- Parameters:
max_num_components – Maximum number of components
weight_concentration_prior (
Optional
[float
]) – The dirichlet concentration of each component on the weight distribution (Dirichlet), default to None, optional This is commonly called gamma in the literature. The higher concentration puts more mass in the center and will lead to more components being active, while a lower concentration parameter will lead to more mass at the edge of the mixture weights simplex. The value of the parameter must be greater than 0. If it is None, it’s set to 1. / n_components.covariance_type (
str
) – {‘full’, ‘tied’, ‘diag’, ‘spherical’}, defaults to “full”, optional String describing the type of covariance parameters to use. Must be one of:'full'
: each component has its own general covariance matrix.'tied'
: all components share the same general covariance matrix.'diag'
: each component has its own diagonal covariance matrix.'spherical'
: each component has its own single variance.n_init (
int
) – number of initialization to perform, defaults to 1, optional The result with the highest lower bound value on the likelihood is kept.tol (
float
) – tolerance threshold, defaults to 1.0e-3, optional EM iterations will stop when the lower bound average gain on the likelihood (of the training data with respect to the model) is below this threshold.max_iter (
int
) – maximum number of EM iterations to perform, defaults to 100, optionalrandom_state (
Optional
[int
]) – random seed given to the method chosen to initialize the parameters, defaults to None, optionalinit_params (
str
) – {‘kmeans’, ‘k-means++’, ‘random’, ‘random_from_data’}, default=’kmeans’, optional The method used to initialize the weights, the means and the covariances. String must be one of:'kmeans'
: responsibilities are initialized using kmeans.k-means++
: responsibilities are initialized using a k-means++ algorithm.'random'
: responsibilities are initialized randomly.'random_from_data'
: responsibilities are initialized randomly from the data.warm_start (
bool
) – flag to use warm start, defaults to False, optional If ‘warm_start’ is True, the solution of the last fitting is used as initialization for the next call of fit(). This can speed up convergence when fit is called several times on similar problems. See the Glossary.expand_factor (
int
) – factor used when converting the ensemble from weighted to unweighted, defaults to 10, optionalslice_sampling (
bool
) – flag to use slice sampling, defaults to False, optional
- draw_samples_within_bounds(system, num=1)
Draw new parameter samples within the user-defined upper and lower bounds
- Parameters:
system (
Type
[DynamicSystem
]) – Dynamic system classnum (
int
) – Number of samples to draw
- Returns:
New parameter samples
- expand_and_normalize_weighted_samples(weights, system)
Converting an ensemble of samples with importance weights into an ensemble of samples with equal weights.
For instance, a sample is duplicated a few times depending on the importance weight associated to it.
- Parameters:
weights (
ndarray
) – Importance weights associated with the ensemblesystem (
Type
[DynamicSystem
]) – Dynamic system class
- Returns:
Expanded unweighted samples
- classmethod from_dict(obj)
Initialize the class using a dictionary style
- Parameters:
obj (
dict
) – a dictionary containing the keys and values to construct a GMM object- Returns:
a GMM object
- load_gmm_from_file(proposal_data_file='proposal_density.pkl')
Load the Gaussian mixture model from a file.
- Parameters:
proposal_data_file (
str
) – Name of the file that stores the trained Gaussian mixture model
- regenerate_params(weight, system)
Regenerating new samples by training and sampling from a Gaussian mixture model.
- Parameters:
weight (
ndarray
) – Posterior found by the data assimilationsystem (
Type
[DynamicSystem
]) – Dynamic system class
- Returns:
Resampled parameter data
- save_gmm_to_file(proposal_data_file='proposal_density.pkl')
Save the Gaussian mixture model to a file.
- Parameters:
proposal_data_file (str) –
- train(weight, system)
Train the Gaussian mixture model.
- Parameters:
weight (
ndarray
) – Posterior found by the data assimilationsystem (
Type
[DynamicSystem
]) – Dynamic system class
-
covariance_type:
str
= 'full' String describing the type of covariance parameters to use.
-
expand_factor:
int
= 10 the factor used when converting and populating the ensemble from weighted to unweighted, defaults to 10.
-
gmm:
Type
[BayesianGaussianMixture
] The class of the Gaussian Mixture Model
-
init_params:
str
= 'kmeans' The method used to initialize the weights, the means and the covariances.
-
max_iter:
int
= 100 maximum number of EM iterations to perform, defaults to 100
-
max_num_components:
int
= 0 Maximum number of components
- max_params = None
Current maximum values of the parameters
-
n_init:
int
= 1 number of initialization to perform, defaults to 1.
-
random_state:
int
random seed given to the method chosen to initialize the weights, the means and the covariances.
- slice_sampling: False
flag to use slice sampling, defaults to False.
-
tol:
float
= 0.001 tolerance threshold, defaults to 1.0e-3.
-
warm_start:
bool
= False flag to use warm start, defaults to False.
-
weight_concentration_prior:
float
= 0.0 The dirichlet concentration of each component on the weight distribution (Dirichlet), default to None.
- generate_params_qmc(system, num_samples, method='halton', seed=None)
This is the class to uniformly draw samples in n-dimensional space from a low-discrepancy sequence or a Latin hypercube.
See Quasi-Monte Carlo.
- Parameters:
system (
Type
[DynamicSystem
]) – Dynamic system classnum_samples (
int
) – Number of samples to drawmethod (
str
) – Method to use for Quasi-Monte Carlo sampling. Options are “halton”, “sobol”, and “LH”seed – Seed for the random number generator
RNN Module
data pre-processing
- class Preprocessor
Bases:
ABC
- classmethod get_dimensions(data)
Extract dimensions of sample from a tensorflow dataset.
- Parameters:
data (
DatasetV2
) – The dataset to extract from.- Returns:
Dictionary containing: sequence_length, num_load_features, num_contact_params, num_labels
- classmethod get_split(dataset, inds)
Given a dataset and the indexes, return the sub-dataset containing only the elements of indexes in
inds
. The returned dataset respects the format ofinputs
andoutputs
(see more info in return).- Parameters:
dataset (
tuple
) – Full dataset to split on, in form of a tuple(inputs, outputs)
, whereinputs
is a dictionary and its keys andoutputs
are numpy arrays.inds (
array
) – Indexes of the elements indataset
that are going to be gathered in this specific split.
- Returns:
tuple of the split, 2 dimensions: * inputs: dictionary containing
'load_sequence'
and'contact_params'
. * outputs: Labels
- classmethod make_splits(dataset, train_frac, val_frac, seed)
Split data into training, validation, and test sets. The split is done on a sample by sample basis, so sequences are not broken up. It is done randomly across all pressures and experiment types present.
- Parameters:
dataset (
tuple
) – Full dataset to split on, in form of a tuple (inputs, outputs), where inputs is a dictionary and its keys and the outputs are numpy arrays.train_frac (
float
) – Fraction of data used for training set.val_frac (
float
) – Fraction of data used for validation set. Test fraction is the remaining.seed (
int
) – Random seed used to make the split.
- Returns:
Dictionary containing
'train'
,'val'
, and'test'
datasets.
- classmethod pad_initial(array, pad_length, axis=1)
Add
pad_length
copies of the initial state in the sequence to the start. This is used to be able to predict also the first timestep from a window of the same size.- Parameters:
array (
array
) – Array that is going to be modifiedpad_lenght – number of copies of the initial state to be added at the beggining of
array
.pad_length (int) –
- Returns:
Modified array
- abstract prepare_datasets(**kwargs)
Abstract class to be implemented in each of the child classes. Must return a tuple containing:
data_split
: Dictionary with keys'train'
,'val'
,'test'
, and values of the corresponding tensorflow Datasets.train_stats
: Dictionary containing the shape of the data:sequence_length
,num_input_features
,num_contact_params
,num_labels
, and'mean'
and'std'
of the training set, in casestandardize_outputs
is True.
- abstract prepare_single_dataset()
Abstract class to be implemented in each of the child classes. Must return a Tensorflow dataset containing a Tuple (inputs, outputs)
inputs
: Dictionary with keys'load_sequence'
,'contact_parameters'
, that are the corresponding tensorflow Datasets to input to an rnn model.outputs
: tensorflow Dataset containing the outputs or labels.
- classmethod standardize_outputs(split_data)
Standardize outputs of
split_data
using the mean and std of the training data taken over both the samples and the timesteps. The 3 datasets'train'
,'val'
,'test'
will be standardized.- Parameters:
split_data (
dict
) – dictionary containing'train'
,'val'
,'test'
keys and the respective datasets.- Returns:
Tuple containing:
standardized_splits
: Dictionary containing the standardized datasets.train_stats
: Dictionary with the metrics (mean and standard deviation) used to standardize the data.
- classmethod warning_config_field(key, config, default, add_default_to_config=False)
Raises a warning if key is not included in config dictionary. Also informs the default value that will be used. If add_default_to_config=True, then it adds the key and its default value to config.
- Parameters:
key (
str
) – string of the key in dictionaryconfig
to look for.config (
dict
) – Dictionary containing keys and valuesdefault – default value associated to the key in config.
add_default_to_config (
bool
) – Wether to add the default value and the key to config in case it does not exists.
- Returns:
Updated
config
dictionary.
- class PreprocessorTriaxialCompression(**kwargs)
Bases:
Preprocessor
Class to Preprocess data of triaxial compression experiments, inheriting from abstract class
Preprocessor
Attributes: see
PreprocessorTriaxialCompression.get_default_config()
- classmethod check_config(config)
Checks that values requiring an input from the user would be specified in config.
- Parameters:
config (
dict
) – Dictionary containing the values of different arguments.- Returns:
Updated config dictionary.
- classmethod get_default_config()
Returns a dictionary with default values for the configuration of data preparation. Possible fields are:
'raw_data'
: Path to hdf5 file generated using parse_data_YADE.py'pressure'
and'experiment_type'
: Name of the subfield of dataset to consider. It can also be ‘All’.'standardize_outputs'
: If True transform the data labels to have zero mean and unit variance. Also, in train_stats the mean and variance of each label will be stored, so that can be used to transform predicitons. (This is very usful if the labels are not between [0,1])'add_e0'
: Whether to add the initial void ratio (output) as a contact parameter.'add_pressure'
: Wether to add the pressure to contact_parameters.'add_experiment_type'
: Wether to add the experiment type to contact_parameters.'train_frac'
: Fraction of the data used for training, between [0,1].'val_frac'
: Fraction of the data used for validation, between [0,1]. The fraction of the data used for test is then1 - train_frac - val_frac
.'window_size'
: int, number of steps composing a window.'window_step'
: int, number of steps between consecutive windows (default = 1).'pad_length'
: int, equals towindow_size
. Length of the sequence that with be pad at the start.
- Returns:
Dictionary containing default values of the arguments that the user can set.
- classmethod get_dimensions(data)
Extract dimensions of sample from a tensorflow dataset.
- Parameters:
data (
DatasetV2
) – The dataset to extract from.- Returns:
Dictionary containing: sequence_length, num_load_features, num_contact_params, num_labels
- classmethod get_split(dataset, inds)
Given a dataset and the indexes, return the sub-dataset containing only the elements of indexes in
inds
. The returned dataset respects the format ofinputs
andoutputs
(see more info in return).- Parameters:
dataset (
tuple
) – Full dataset to split on, in form of a tuple(inputs, outputs)
, whereinputs
is a dictionary and its keys andoutputs
are numpy arrays.inds (
array
) – Indexes of the elements indataset
that are going to be gathered in this specific split.
- Returns:
tuple of the split, 2 dimensions: * inputs: dictionary containing
'load_sequence'
and'contact_params'
. * outputs: Labels
- classmethod make_splits(dataset, train_frac, val_frac, seed)
Split data into training, validation, and test sets. The split is done on a sample by sample basis, so sequences are not broken up. It is done randomly across all pressures and experiment types present.
- Parameters:
dataset (
tuple
) – Full dataset to split on, in form of a tuple (inputs, outputs), where inputs is a dictionary and its keys and the outputs are numpy arrays.train_frac (
float
) – Fraction of data used for training set.val_frac (
float
) – Fraction of data used for validation set. Test fraction is the remaining.seed (
int
) – Random seed used to make the split.
- Returns:
Dictionary containing
'train'
,'val'
, and'test'
datasets.
- classmethod pad_initial(array, pad_length, axis=1)
Add
pad_length
copies of the initial state in the sequence to the start. This is used to be able to predict also the first timestep from a window of the same size.- Parameters:
array (
array
) – Array that is going to be modifiedpad_lenght – number of copies of the initial state to be added at the beggining of
array
.pad_length (int) –
- Returns:
Modified array
- prepare_datasets(seed=42)
Convert raw data into preprocessed split datasets. First split the data into train, val and test datasets and then apply the Sliding windows transformation. This is to avoid having some parts of a dataset in train and some in val and/or in test (i.e. data leak).
- Parameters:
seed (
int
) – Random seed used to split the datasets.- Returns:
Tuple (split_data, train_stats)
split_data
: Dictionary with keys'train'
,'val'
,'test'
, and values the corresponding tensorflow Datasets.train_stats
: Dictionary containing the shape of the data:sequence_length
,num_load_features
,num_contact_params
,num_labels
, and'mean'
and'std'
of the training set, in casestandardize_outputs
is True.
- prepare_single_dataset()
Convert raw data into a tensorflow dataset with compatible format to predict and evaluate a rnn model.
- Returns:
Tensorflow dataset containing a Tuple (inputs, outputs)
inputs
: Dictionary with keys'load_sequence'
,'contact_parameters'
, that are the corresponding tensorflow Datasets to input to an rnn model.outputs
: tensorflow Dataset containing the outputs or labels.
- classmethod standardize_outputs(split_data)
Standardize outputs of
split_data
using the mean and std of the training data taken over both the samples and the timesteps. The 3 datasets'train'
,'val'
,'test'
will be standardized.- Parameters:
split_data (
dict
) – dictionary containing'train'
,'val'
,'test'
keys and the respective datasets.- Returns:
Tuple containing:
standardized_splits
: Dictionary containing the standardized datasets.train_stats
: Dictionary with the metrics (mean and standard deviation) used to standardize the data.
- classmethod warning_config_field(key, config, default, add_default_to_config=False)
Raises a warning if key is not included in config dictionary. Also informs the default value that will be used. If add_default_to_config=True, then it adds the key and its default value to config.
- Parameters:
key (
str
) – string of the key in dictionaryconfig
to look for.config (
dict
) – Dictionary containing keys and valuesdefault – default value associated to the key in config.
add_default_to_config (
bool
) – Wether to add the default value and the key to config in case it does not exists.
- Returns:
Updated
config
dictionary.
models
Module containing a function that creates a RNN model.
- rnn_model(input_shapes, window_size=20, lstm_units=50, dense_units=20, seed=42, **_)
Neural network with an LSTM layer.
Takes in a load sequence and contact parameters, and outputs the macroscopic responses. The contact parameters are used to initialize the hidden state of the LSTM.
- Parameters:
input_shapes (
dict
) – Dictionary containing ‘num_load_features’, ‘num_contact_params’, ‘num_labels’. It can contain other keys but hese are the ones used here.window_size (
int
) – Length of time window.lstm_units (
int
) – Number of units of the hidden state of the LSTM.dense_units (
int
) – Number of units used in the dense layer after the LSTM.seed (
int
) – The random seed used to initialize the weights.
- Returns:
A Keras model.
windows
- extract_tensors(data)
Given a tensorflow Dataset extract all tensors.
- Parameters:
data (
DatasetV2
) – Tensorflow dataset.- Returns:
3 numpy arrays: inputs, contacts, outputs.
- predict_over_windows(inputs, model, window_size, sequence_length)
Take a batch of full sequences, iterate over windows making predictions. It splits up the sequence into windows of given length, each offset by one timestep, uses the model to make predictions on all of those windows, and concatenates the result into a whole sequence again. Note the length of the output sequence will be shorter by the window_size than the input sequence.
- Parameters:
inputs (
dict
) – dict containing inputs: ‘load_sequence’ and ‘contact_parameters’, both being tensorflow.Tensor.model (
Model
) – The model to predict with.window_size (
int
) – Number of timesteps in a single window.sequence_length (
int
) – Number of timesteps in a full sequence.
- Returns:
tensorflow.Tensor of predicted sequences.
- windowize_single_dataset(data, window_size=1, window_step=1, seed=42, **_)
Take a dataset ((inputs, contact_params),outputs) of sequences of shape N, S, L and output another dataset of shorter sequences of size
window_size
, taken at intervalswindow_step
. The resulting dataset will haveM=((sequence_length - window_size)/window_step) + 1)
samples, corresponding to independent windows in the given dataset. For a given window taken from inputs, the output is the element of outputs sequence at the last position of the window.Note
For a clear picture of what goes on here check Sliding windows section in the documentation.
Data is shuffled.
- Parameters:
data (
DatasetV2
) – dataset of sequences of shape N, S, L.window_size (
int
) – Size of the window.window_step (
int
) – Offset between subsequent windows.seed (
int
) – Random seed.
- Returns:
inputs
of shape: (M, window_size,num_load_features
), with M >> N.outputs
of shape: (M, L_outputs)
- windowize_train_val_test(split_data, window_size, window_step, **_)
Convert sequences into windows of given length. Leave test set untouched.
- Parameters:
split_data (
dict
) – Dictionary with keys ‘train’, ‘val’, ‘test’ pointing to tensorflow datasets.window_size (
int
) – Number of timesteps to include in a window.window_step (
int
) – Offset between subsequent windows.
- Returns:
windows: Dictionary of dataset splits, where the training and validation splits have been modified into windows.
train
Train a model to predict macroscopic features of a DEM simulation.
- get_default_config()
Returns a dictionary with default values for the configuration of RNN model and training procedure. Possible fields are:
RNN model
'window_size'
: int, number of steps composing a window.'window_step'
: int, number of steps between consecutive windows (default = 1).'pad_length'
: int, equals towindow_size
. Length of the sequence that with be pad at the start.'lstm_units'
: int, number of neurons or units in LSTM layer.'dense_units'
: int, number of neurons or units of dense layer.
Training procedure
'patience'
: patience of tf.keras.callbacks.EarlyStopping.'epochs'
: Maximum number of epochs.'learning_rate'
: double, learning_rate of tf.keras.optimizers.Adam.'batch_size'
: Size of the data batches per training step.'save_weights_only'
: BooleanTrue: Only the weights will be saved (Recommended fro compatibility across platforms).
False: The whole model will be saved.
- Returns:
Dictionary containing default values of the arguments that the user can set.
- train(preprocessor, config=None, model=None)
Train a model and report to weights and biases.
If called in the framework of a sweep: A sweep can be created from the command line using a configuration file, for example example_sweep.yaml, as:
wandb sweep example_sweep.yaml
And run with the line shown subsequently in the terminal. The config is loaded from the yaml file.- Parameters:
preprocessor (
Preprocessor
) – Preprocessor object to load and prepare the data.config – dictionary containing model and training configurations.
model (
Optional
[Model
]) – Keras model ifNone
is passed then anrnn_model
will be created (default).
- Returns:
Same as
tf.keras.Model.fit()
: A History object. Its History.history attribute is a record of training loss values and metrics values at successive epochs, as well as validation loss values and validation metrics values.
- train_without_wandb(preprocessor, config=None, model=None)
Train a model locally: no report to wandb. Saves either the model or its weight to folder outputs.
- Parameters:
preprocessor (
Preprocessor
) – Preprocessor object to load and prepare the data.config – dictionary containing taining hyperparameters and some model parameters.
model (
Optional
[Model
]) – Keras model ifNone
is passed then anrnn_model
will be created (default).
- Returns:
Same as
tf.keras.Model.fit()
: A History object. Its History.history attribute is a record of training loss values and metrics values at successive epochs, as well as validation loss values and validation metrics values.
predict
Module containing functions to load a trained RNN model and make a prediction.
- get_best_run_from_sweep(entity_project_sweep_id, preprocessor=None)
Load the best performing model found with a weights and biases sweep. Also load the data splits it was trained on.
- Parameters:
entity_project_sweep_id (
str
) – string of form <user>/<project>/<sweep_id>preprocessor (
Optional
[Preprocessor
]) – Preprocessor object to load and prepare the data. If None is given, then a PreprocessorTriaxialCompression will be considered
- Returns:
model: The trained model with the lowest validation loss.
data: The train, val, test splits as used during training.
stats: Some statistics on the training set and configuration.
config: The configuration dictionary used to train the model.
- get_pretrained_model(path_to_model)
Loads configuration, training statistics and model of a pretrained model.
Reads train_stats, and creates dataset.
- Parameters:
path_to_model (
str
) – str or pathlib.Path to the folder where is stored.- Returns:
model: keras model ready to use.
train_stats: Array containing the values used to standardize the data (if config.standardize_outputs = True), and lenghts of sequences, load_features, contact_params, labels, window_size and window_step.
config: dictionary with the model configuration
- load_config(path_to_model)
Searches for the configuration (of the model training) file in ‘path_to_model’. Read config.yaml into a python dictionary equivalent to config. config.yaml contains information about hyperparameters and model parameters, is generated in every run of wandb. Raises FileNotFoundError if there are not files matching possible formats.
- Parameters:
path_to_model (
Path
) – str or pathlib.Path to the folder where is stored.- Returns:
Dictionary config.
- load_model(path_to_model, train_stats, config)
Searches for the file containing the saved model in ‘path_to_model’. Raises FileNotFoundError if there are not files matching possible formats.
- Parameters:
path_to_model (
Path
) – str or pathlib.Path to the folder where is stored.train_stats (
dict
) – Dictionary containing different dimensions of the data used to trained and standardization values. Its contents are saved in ‘path_to_model/train_stats.npy’.config (
dict
) – Dictionary containing the configuration of the training for such model.
- Returns:
Keras model.
- predict_macroscopics(model, data, train_stats, config, batch_size=256)
Use the given model to predict the features of the given data. If ‘standardize_outputs’ in config, rescale the predictions to their original units.
- Parameters:
model (
Model
) – Keras RNN modeldata (
DatasetV2
) – Tensorflow dataset containing inputs: ‘load_sequence’ and ‘contact_parameters’, and outputs.train_stats (
dict
) – Dictionary containing statistics of the training set.config (
dict
) – Dictionary containing the configuration with which the model was trained.batch_size (
int
) – Size of batches to use.
- Returns:
predictions: tf.Tensor containing the predictions in original units.