Dakota
Version 6.20
Explore and Predict with Confidence
|
Public Member Functions | |
ExperimentData () | |
default constructor | |
ExperimentData (const ProblemDescDB &prob_desc_db, const SharedResponseData &srd, short output_level) | |
typical DB-based constructor | |
ExperimentData (size_t num_experiments, size_t num_config_vars, const boost::filesystem::path &data_prefix, const SharedResponseData &srd, const StringArray &variance_types, short output_level, std::string scalarDataFilename="") | |
temporary? constructor for testing | |
ExperimentData (size_t num_experiments, const SharedVariablesData &svd, const SharedResponseData &srd, const VariablesArray &configVars, const IntResponseMap &all_responses, short output_level) | |
Bayesian experimental design constructor. Passed SVD has calibration parameters as active and config vars as inactive variables. Passed configVars have config vars as active. More... | |
void | load_data (const std::string &context_message, const Variables &vars_with_state_as_config) |
Load experiments from data files (simple scalar or field) | |
void | add_data (const SharedVariablesData &svd, const Variables &one_configvars, const Response &one_response) |
Add one data point to the experimental data set. Used for Bayesian experimental design. Passed SVD has calibration parameters as active and config vars as inactive variables. Passed configVars have config vars as active. | |
size_t | num_experiments () const |
retrieve the number of experiments | |
size_t | num_total_exppoints () const |
retrieve the total number of experimental data points over all experiments | |
size_t | num_scalar_primary () const |
retrieve the number of scalars (applies to all experiments) | |
size_t | num_fields () const |
retrieve the number of fields (applies to all experiments) | |
size_t | num_config_vars () const |
number of configuration variables | |
std::vector< RealVector > | config_vars_as_real () const |
values of the configuration variables, 1 RealVector per experiment omits string variables as historically used in NonDBayes More... | |
const std::vector< Variables > & | configuration_variables () const |
std::vector< Variables > & | configuration_variables () |
const RealVector & | all_data (size_t experiment) |
return contiguous vector of all data (scalar, followed by field) for the specified experiment | |
const Response & | response (size_t experiment) |
return response for the specified experiment | |
void | per_exp_length (IntVector &per_length) const |
return the individual sizes of the experimental data lengths (all function values, scalar and field) | |
const IntVector & | field_lengths (size_t experiment) const |
return the field lengths for specified experiment index | |
Real | scalar_data (size_t response, size_t experiment) |
retrieve the data value for the given response, for the given experiment | |
RealVector | field_data_view (size_t response, size_t experiment) const |
retrieve a view of the field data for the given response, for the given experiment | |
RealMatrix | field_coords_view (size_t response, size_t experiment) const |
retrieve a view of the field data coordinates for the given response, for the given experiment | |
void | fill_primary_function_labels (StringArray &expanded_labels) const |
populate the passed array with num_total_exppoints labels | |
bool | variance_type_active (short variance_type) const |
whether the specified variance type (enum value) is present and active | |
bool | variance_active () const |
whether any variance type is active | |
Real | apply_covariance (const RealVector &residuals, size_t experiment) const |
apply the covariance responses to compute the triple product v'*inv(C)*v for the given experiment | |
void | apply_covariance_inv_sqrt (const RealVector &residuals, size_t experiment, RealVector &weighted_residuals) const |
apply inverse sqrt of the covariance to compute weighted residuals | |
void | apply_covariance_inv_sqrt (const RealMatrix &gradients, size_t experiment, RealMatrix &weighted_gradients) const |
apply inverse sqrt of the covariance to compute weighted gradients | |
void | apply_covariance_inv_sqrt (const RealSymMatrixArray &hessians, size_t experiment, RealSymMatrixArray &weighted_hessians) const |
apply inverse sqrt of the covariance to compute weighted Hessians | |
void | apply_simulation_error (const RealVector &simulation_error, size_t experiment) |
apply simulation error to experiment data | |
void | get_main_diagonal (RealVector &diagonal, size_t experiment) const |
return a (copy) vector containing the main diagonal entries of a specified experimental covariance matrix | |
void | cov_std_deviation (RealVectorArray &std_deviation) const |
get the standard deviation of the observation error process, one vector per experiment | |
void | cov_as_correlation (RealSymMatrixArray &corr_matrix) const |
get the observation error covariance as a correlation matrix, one vector per experiment | |
void | covariance (int exp_ind, RealSymMatrix &cov_mat) const |
retrieve an individual covariance entry as a dense matrix | |
void | form_residuals (const Response &sim_resp, Response &residual_resp) const |
form residuals for all experiments, interpolating if necessary; one simulation response maps to all experiments More... | |
void | form_residuals (const Response &sim_resp, const size_t curr_exp, Response &residual_resp) const |
Populate the portion of residual_resp corresponding to experiment curr_exp; the passed simulation response maps only to the specified experiment. | |
void | form_residuals (const Response &sim_resp, size_t exp_num, const ShortArray &total_asv, size_t residual_resp_offset, Response &residual_resp) const |
form residuals for an individual experiment, interpolating if necessary More... | |
void | recover_model (size_t num_pri_fns, RealVector &model_fns) const |
recover original model from the first experiment block in a full set of residuals; works in no interpolation case only (sizes same) More... | |
bool | interpolate_flag () const |
flag for interpolation. If 0, no interpolation. If 1, interpolate. | |
void | interpolate_simulation_data (const Response &sim_resp, size_t exp_num, const ShortArray &total_asv, size_t exp_offset, Response &interp_resp) const |
Interpolate simulation data (values, gradients and hessians) onto the coordinates of the experimental data. | |
void | scale_residuals (const Response &residual_response, RealVector &scaled_residuals) const |
Apply the experiment data covariance to the residual data (scale functions by Gamma_d^{-1/2}), returning in scaled_residuals. | |
void | scale_residuals (Response &residual_response) const |
Apply the experiment data covariance to the residual data in-place (scale functions, gradients, and Hessians by Gamma_d^{-1/2}) | |
void | build_gradient_of_sum_square_residuals (const Response &resp, RealVector &ssr_gradient) |
Build the gradient of the ssr from residuals and function gradients based on the response's active set request vector. | |
void | build_gradient_of_sum_square_residuals (const Response &resp, const ShortArray &asrv, RealVector &ssr_gradient) |
Build the gradient of the ssr from residuals and function gradients using the passed active set request vector (overrides the response's request vector) | |
void | build_gradient_of_sum_square_residuals_from_response (const Response &resp, const ShortArray &asrv, int exp_ind, RealVector &ssr_gradient) |
Update the gradient of ssr with the values from the gradient associated with a single experiment. | |
void | build_gradient_of_sum_square_residuals_from_function_data (const RealMatrix &func_gradients, const RealVector &residuals, RealVector &ssr_gradient, const ShortArray &asrv) |
Construct the gradient of the sum of squares of residuals. More... | |
void | build_hessian_of_sum_square_residuals (const Response &resp, RealSymMatrix &ssr_hessian) |
Build the hessian of the ssr from residuals, function gradients and function hessians based on the response's active set request vector. | |
void | build_hessian_of_sum_square_residuals (const Response &resp, const ShortArray &asrv, RealSymMatrix &ssr_hessian) |
Build the hessian of the ssr from residuals, function gradients and function hessians using the passed active set request vector (overrides the response's request vector) | |
void | build_hessian_of_sum_square_residuals_from_response (const Response &resp, const ShortArray &asrv, int exp_ind, RealSymMatrix &ssr_hessian) |
Update the hessian of ssr with the values from the hessian associated with a single experiment. | |
void | build_hessian_of_sum_square_residuals_from_function_data (const RealSymMatrixArray &func_hessians, const RealMatrix &func_gradients, const RealVector &residuals, RealSymMatrix &ssr_hessian, const ShortArray &asrv) |
Construct the hessian of the sum of squares of residuals. More... | |
void | scale_residuals (const RealVector &multipliers, unsigned short multiplier_mode, size_t num_calib_params, Response &residual_response) const |
in-place scale the residual response (functions, gradients, Hessians) by sqrt(multipliers), according to blocks indicated by multiplier mode More... | |
Real | cov_determinant (const RealVector &multipliers, unsigned short multiplier_mode) const |
returns the determinant of (covariance block-scaled by the passed multipliers) More... | |
Real | half_log_cov_determinant (const RealVector &multipliers, unsigned short multiplier_mode) const |
returns the log of the determinant of (covariance block-scaled by the passed multipliers) More... | |
void | half_log_cov_det_gradient (const RealVector &multipliers, unsigned short multiplier_mode, size_t hyper_offset, RealVector &gradient) const |
populated the passed gradient with derivatives w.r.t. the hyper-parameter multipliers, starting at hyper_offset (must be sized) More... | |
void | half_log_cov_det_hessian (const RealVector &multipliers, unsigned short multiplier_mode, size_t hyper_offset, RealSymMatrix &hessian) const |
populated the passed Hessian with derivatives w.r.t. the hyper-parameter multipliers, starting at hyper_offset (must be sized) More... | |
StringArray | hyperparam_labels (unsigned short multiplier_mode) const |
generate variable labels for the covariance (error) multiplier hyperparams | |
Protected Member Functions | |
ShortArray | determine_active_request (const Response &resid_resp) const |
Perform check on the active request vector to make sure it is amenable to interpolation of simulation data and application of apply covariance. | |
SizetArray | residuals_per_multiplier (unsigned short multiplier_mode) const |
count the number of residuals influenced by each multiplier More... | |
void | generate_multipliers (const RealVector &multipliers, unsigned short multiplier_mode, RealVector &expanded_multipliers) const |
Generate a set of multipliers commensurate with the residual size for the total experiment data set. Instead of repeating the loops all over the place, generate an expanded set of multipliers; the conditionals get too complicated otherwise. | |
void | resid2mult_map (unsigned short multiplier_mode, IntVector &resid2mult_indices) const |
return the index of the multiplier that affects each residual | |
Private Member Functions | |
void | initialize (const StringArray &variance_types, const SharedResponseData &srd) |
shared body of constructor initialization | |
void | parse_sigma_types (const StringArray &sigma_types) |
parse user-provided sigma type strings and populate enums More... | |
void | update_data_properties () |
After constructing or adding data, update properties like experiment lengths, determinant, etc. | |
void | load_experiment (size_t exp_index, std::ifstream &scalar_data_stream, size_t num_field_sigma_matrices, size_t num_field_sigma_diagonals, size_t num_field_sigma_scalars, size_t num_field_sigma_none, Response &exp_resp) |
Load a single experiment exp_index into exp_resp. More... | |
void | read_scalar_sigma (std::ifstream &scalar_data_stream, RealVector &sigma_scalars, IntVector &scalar_map_indices) |
read or default populate the scalar sigma | |
RealVector | residuals_view (const RealVector &residuals, size_t experiment) const |
Return a view (to allowing updaing in place) of the residuals associated with a given experiment, from a vector contaning residuals from all experiments. | |
RealMatrix | gradients_view (const RealMatrix &gradients, size_t experiment) const |
Return a view (to allowing updaing in place) of the gradients associated with a given experiment, from a matrix contaning gradients from all experiments. | |
RealSymMatrixArray | hessians_view (const RealSymMatrixArray &hessians, size_t experiment) const |
Return a view (to allowing updaing in place) of the hessians associated with a given experiment, from an array contaning the hessians from all experiments. | |
Private Attributes | |
bool | calibrationDataFlag |
whether the user specified a calibration data block | |
size_t | numExperiments |
the total number of experiments | |
size_t | numConfigVars |
number of configuration (state) variables to read for each experiment | |
UShortArray | varianceTypes |
type of variance specified for each variable, one per response group; empty varianceType indicates none specified by user | |
Real | covarianceDeterminant |
cached product of each experiment covariance's determinant | |
Real | logCovarianceDeterminant |
cached sum of each experiment covariance's log determinant | |
boost::filesystem::path | dataPathPrefix |
path to prepend to any data file names | |
String | scalarDataFilename |
the user-specied scalar data filename | |
unsigned short | scalarDataFormat |
tabular format of the simple scalar data file; supports TABULAR_NONE, TABULAR_HEADER, TABULAR_EVAL_ID, TABULAR_EXPER_ANNOT | |
size_t | scalarSigmaPerRow |
number of sigma values to read from each row in simple data file format (calculated from variance types strings | |
bool | readSimFieldCoords |
whether to read coordinate data files for simulation fields | |
SharedResponseData | simulationSRD |
archived shared data for use in sizing fields, total functions (historically we read all functions, including constraints, which might not be correct) | |
bool | interpolateFlag |
flag for interpolation. | |
short | outputLevel |
output verbosity level | |
std::vector< Response > | allExperiments |
Vector of numExperiments ExperimentResponses, holding the observed data and error (sigma/covariance) for each experiment. | |
std::vector< Variables > | allConfigVars |
Vector of numExperiments configurations at which data were gathered; empty if no configurations specified. The inactive state variables are used to store the configuration settings. | |
IntVector | experimentLengths |
Length of each experiment. | |
IntVector | expOffsets |
function index offsets for individual experiment data sets | |
The ExperimentData class is used to read and populate data (currently from user-specified files and/or the input spec) relating to experimental (physical observations) data for the purposes of calibration. Such data may include (for example): number of experiments, configuration variables, type of data (scalar vs. functional), treatment of sigma (experimental uncertainties). This class also provides an interpolation capability to interpolate between simulation or experimental data so that the differencing between simulation and experimental data may be performed properly.
ExperimentData | ( | size_t | num_experiments, |
const SharedVariablesData & | svd, | ||
const SharedResponseData & | srd, | ||
const VariablesArray & | config_vars, | ||
const IntResponseMap & | all_responses, | ||
short | output_level | ||
) |
Bayesian experimental design constructor. Passed SVD has calibration parameters as active and config vars as inactive variables. Passed configVars have config vars as active.
Used in Hi2Lo Bayesian experimental design; passed config vars are active, but stored here as inactive.
References ExperimentData::allConfigVars, ExperimentData::allExperiments, SharedResponseData::copy(), Response::copy(), ExperimentData::numConfigVars, ExperimentData::numExperiments, ExperimentData::outputLevel, SharedResponseData::response_type(), ExperimentData::simulationSRD, Dakota::size_and_fill(), Dakota::svd(), Response::update(), and ExperimentData::update_data_properties().
std::vector< RealVector > config_vars_as_real | ( | ) | const |
values of the configuration variables, 1 RealVector per experiment omits string variables as historically used in NonDBayes
Skips string vars rather than converting to indices
References ExperimentData::allConfigVars, Dakota::copy_data_partial(), and Dakota::merge_data_partial().
Referenced by NonDGPMSABayesCalibration::fill_experiment_data().
form residuals for all experiments, interpolating if necessary; one simulation response maps to all experiments
This assumes the souce gradient/Hessian are size less or equal to the destination response, and that the leading part is to be populated.
References ExperimentData::determine_active_request(), ExperimentData::numExperiments, and ExperimentData::per_exp_length().
Referenced by DataTransformModel::archive_submodel_responses(), DataTransformModel::derived_evaluate(), ExperimentData::form_residuals(), DataTransformModel::primary_resp_differencer(), and DataTransformModel::transform_response_map().
void form_residuals | ( | const Response & | sim_resp, |
size_t | exp_ind, | ||
const ShortArray & | total_asv, | ||
size_t | exp_offset, | ||
Response & | residual_resp | ||
) | const |
form residuals for an individual experiment, interpolating if necessary
This assumes the souce gradient/Hessian are size less or equal to the destination response, and that the leading part is to be populated.
References ExperimentData::allExperiments, ExperimentData::field_data_view(), Response::field_lengths(), Response::function_gradient_view(), Response::function_gradients(), Response::function_gradients_view(), Response::function_hessian_view(), Response::function_hessians(), Response::function_hessians_view(), Response::function_values(), Response::function_values_view(), ExperimentData::gradients_view(), ExperimentData::hessians_view(), ExperimentData::interpolate_simulation_data(), ExperimentData::interpolateFlag, ExperimentData::num_fields(), ExperimentData::num_scalar_primary(), and ExperimentData::outputLevel.
void recover_model | ( | size_t | num_pri_fns, |
RealVector & | best_fns | ||
) | const |
recover original model from the first experiment block in a full set of residuals; works in no interpolation case only (sizes same)
Add the data back to the residual to recover the model, for use in surrogated-based LSQ where DB lookup will fail (need approx eval DB). best_fns contains primary and secondary responses
References Dakota::abort_handler(), ExperimentData::allExperiments, Response::function_value(), ExperimentData::interpolateFlag, SharedResponseData::num_primary_functions(), and Response::shared_data().
void build_gradient_of_sum_square_residuals_from_function_data | ( | const RealMatrix & | func_gradients, |
const RealVector & | residuals, | ||
RealVector & | ssr_gradient, | ||
const ShortArray & | asrv | ||
) |
Construct the gradient of the sum of squares of residuals.
func_gradients | A matrix containing the gradients of the residual vector |
residuals | A vector of residuals (mismatch between experimental data and the corresponding function values |
asrv | The active set request vector |
Referenced by ExperimentData::build_gradient_of_sum_square_residuals_from_response().
void build_hessian_of_sum_square_residuals_from_function_data | ( | const RealSymMatrixArray & | func_hessians, |
const RealMatrix & | func_gradients, | ||
const RealVector & | residuals, | ||
RealSymMatrix & | ssr_hessian, | ||
const ShortArray & | asrv | ||
) |
Construct the hessian of the sum of squares of residuals.
func_hessians | A list of matrices containing the Hessians of the function elements in the residual vector |
func_gradients | A matrix containing the gradients of the residual vector |
residuals | A vector of residuals (mismatch between experimental data and the corresponding function values |
asrv | The active set request vector |
Referenced by ExperimentData::build_hessian_of_sum_square_residuals_from_response().
void scale_residuals | ( | const RealVector & | multipliers, |
unsigned short | multiplier_mode, | ||
size_t | num_calib_params, | ||
Response & | residual_response | ||
) | const |
in-place scale the residual response (functions, gradients, Hessians) by sqrt(multipliers), according to blocks indicated by multiplier mode
In-place scaling of residual response by hyper-parameter multipliers
References Dakota::abort_handler(), Response::active_set_request_vector(), Response::function_gradient_view(), Response::function_hessian_view(), Response::function_value(), Response::function_value_view(), ExperimentData::num_total_exppoints(), and ExperimentData::resid2mult_map().
Real cov_determinant | ( | const RealVector & | multipliers, |
unsigned short | multiplier_mode | ||
) | const |
returns the determinant of (covariance block-scaled by the passed multipliers)
Determinant of the total covariance used in inference, which has blocks mult_i * I * Cov_i.
References Dakota::abort_handler(), ExperimentData::covarianceDeterminant, ExperimentData::generate_multipliers(), and ExperimentData::num_total_exppoints().
Real half_log_cov_determinant | ( | const RealVector & | multipliers, |
unsigned short | multiplier_mode | ||
) | const |
returns the log of the determinant of (covariance block-scaled by the passed multipliers)
Determinant of half the log of total covariance used in inference, which has blocks mult_i * I * Cov_i.
References Dakota::abort_handler(), ExperimentData::generate_multipliers(), ExperimentData::logCovarianceDeterminant, and ExperimentData::num_total_exppoints().
Referenced by NonDBayesCalibration::log_likelihood(), NonDMUQBayesCalibration::print_results(), and NonDQUESOBayesCalibration::print_results().
void half_log_cov_det_gradient | ( | const RealVector & | multipliers, |
unsigned short | multiplier_mode, | ||
size_t | hyper_offset, | ||
RealVector & | gradient | ||
) | const |
populated the passed gradient with derivatives w.r.t. the hyper-parameter multipliers, starting at hyper_offset (must be sized)
Compute the gradient of scalar f(m) 0.5*log(det(mult*Cov)) w.r.t. mults. Since this is the only use case, we include the 0.5 factor and perform an update in-place.
References ExperimentData::num_total_exppoints(), and ExperimentData::residuals_per_multiplier().
Referenced by NonDBayesCalibration::neg_log_post_resp_mapping().
void half_log_cov_det_hessian | ( | const RealVector & | multipliers, |
unsigned short | multiplier_mode, | ||
size_t | hyper_offset, | ||
RealSymMatrix & | hessian | ||
) | const |
populated the passed Hessian with derivatives w.r.t. the hyper-parameter multipliers, starting at hyper_offset (must be sized)
Compute the gradient of scalar f(m) log(det(mult*Cov)) w.r.t. mults
References ExperimentData::num_total_exppoints(), and ExperimentData::residuals_per_multiplier().
Referenced by NonDBayesCalibration::calculate_evidence(), and NonDBayesCalibration::neg_log_post_resp_mapping().
|
protected |
count the number of residuals influenced by each multiplier
Calculate how many residuals each multiplier affects
References ExperimentData::allExperiments, SharedResponseData::num_field_response_groups(), ExperimentData::num_fields(), SharedResponseData::num_response_groups(), SharedResponseData::num_scalar_primary(), ExperimentData::numExperiments, and ExperimentData::simulationSRD.
Referenced by ExperimentData::half_log_cov_det_gradient(), and ExperimentData::half_log_cov_det_hessian().
|
private |
parse user-provided sigma type strings and populate enums
Validate user-provided sigma specifcation. User can specify 0, 1, or num_response_groups sigmas. If specified, sigma types must be the same for all scalar responses.
References Dakota::abort_handler(), SharedResponseData::num_field_response_groups(), SharedResponseData::num_scalar_primary(), ExperimentData::scalarDataFilename, ExperimentData::scalarSigmaPerRow, ExperimentData::simulationSRD, and ExperimentData::varianceTypes.
Referenced by ExperimentData::initialize().
|
private |
Load a single experiment exp_index into exp_resp.
Load an experiment from a mixture of legacy format data and field data format files
References Dakota::abort_handler(), ExperimentData::dataPathPrefix, Response::field_coords(), Response::field_group_labels(), Response::field_lengths(), ExperimentData::field_lengths(), Response::field_values(), Response::function_labels(), Response::function_value(), Dakota::is_matrix_symmetric(), SharedResponseData::num_field_response_groups(), ExperimentData::num_fields(), SharedResponseData::num_scalar_primary(), ExperimentData::numExperiments, ExperimentData::read_scalar_sigma(), ExperimentData::scalarDataFilename, ExperimentData::scalarSigmaPerRow, Response::set_full_covariance(), ExperimentData::simulationSRD, and ExperimentData::varianceTypes.
Referenced by ExperimentData::load_data().