From 63997b7a1087d3a7a12eea34f4cdb81766b58269 Mon Sep 17 00:00:00 2001 From: Aditya Savara <39929571+AdityaSavara@users.noreply.github.com> Date: Wed, 20 Jan 2021 13:18:25 -0500 Subject: [PATCH] Adding in more Bayesian_dict arguments. (#6) * Update README.rst Moved image assets and added a GUI screenshot * Update README.rst * Moved to Assets Branch * Update fit_fcn.py Replace all for following: change calculate_residuals to calculate_objective_function change calc_resid_output to calc_objective_function_output * Update fit_fcn.py * Update fit_fcn.py changing variable names in fit_fcn and separating the bayesian case into an if statement. * Update fit_fcn.py * Update fit_fcn.py * Update fit_fcn.py Making single return for verbose versus normal case. * Made CheKiPEUQ_from_Frhodo helper module and also put code into fit_fcn Right now, these codes are not working. It is just the beginning. * Update fit_fcn.py * Create ExampleConfig.ini * add comments to fit_fcn.py about what is fed in by "args_list" * Update fit_fcn.py * renaming obs to obs_sim in fit_fcn * Update fit_fcn.py * get_varying_rate_vals_and_bnds * get_varying_rate_vals_and_bnds * pars_uncertainty_distribution code * newOptimization field for creating PE_object. * Added "Force Bayesian" * minor changes * Update CheKiPEUQ_from_Frhodo.py * Update CheKiPEUQ_from_Frhodo.py * Adjusting observed_data and responses shapes for CheKiPEUQ * got the weightings multiplication to work after transposing and transposing back * typo correction * update get_last_obs_sim_interp * moved CheKiPEUQ_PE_object creation into time_adjust_func * Update fit_fcn.py * working on get_log_posterior_density and getting varying_rate_vals * Update fit_fcn.py * forcing 'final' after minimization to be 'residual' * switch to negative log P for objective_function_value * Trying to allow Bayesian way to get past the "QQ" error by feeding residual metrics. * trying 10**neg_logP in lieu of loss_scalar. * adding Bayesian_dict to make code easier to follow * Separating Bayesian into 5 steps to simplify things for Travis * Update fit_fcn.py * Update CheKiPEUQ_integration_notes.txt * Update fit_fcn.py * Update CheKiPEUQ_integration_notes.txt * Create CiteSoftLocal.py * Update CiteSoft call for CheKiPEUQ * Disabling CiteSoft exportations for now * Update CheKiPEUQ_integration_notes.txt * update comment in fit_fcn.py * manually updating fit_fcn.py from Travis's CheKiPEUQ_integration branch Github was not allowing the merge properly, so doing this change manually before merging. * Update CheKiPEUQ_from_Frhodo from the CheKiPEUQ_Integration branch * Update fit_fcn.py * Update fit_fcn.py * creating get_last_obs_sim_interp and also working on Bayesian_dict * Update fit_fcn.py * Update fit_fcn.py Co-authored-by: Travis Sikes <50559900+tsikes@users.noreply.github.com> --- src/optimize/fit_fcn.py | 78 ++++++++++++++++++++++++++++++++--------- 1 file changed, 61 insertions(+), 17 deletions(-) diff --git a/src/optimize/fit_fcn.py b/src/optimize/fit_fcn.py index a0c44b7..d6371c0 100644 --- a/src/optimize/fit_fcn.py +++ b/src/optimize/fit_fcn.py @@ -285,13 +285,21 @@ def __init__(self, input_dict): self.__abort = False if self.opt_settings['obj_fcn_type'] == 'Bayesian': # initialize Bayesian_dictionary if Bayesian selected + #Step 1 of Bayesian: Prepare any variables that need to be passed in for Bayesian PE_object creation. + #Step 2 of Bayesian: populate Bayesian_dict with any variables and uncertainties needed. self.Bayesian_dict = {} + self.Bayesian_dict['pars_uncertainty_distribution'] = self.opt_settings['bayes_dist_type'] #options can be 'Auto', 'Gaussian' or 'Uniform'. # T. Sikes note: optimization is performed on log(rate coefficients), not sure if this is ok or should be np.exp(pars) - self.Bayesian_dict['pars_initial_guess'] = self.x0 + # A. Savara agrees that a log scale is appropriate for rate constants, and suggests base 10, but base e is also fine. + #A. Savara recommends 'uniform' for rate constants and 'gaussian' for things like "log(A)" and "Ea" - self.Bayesian_dict['pars_uncertainty_distribution'] = self.opt_settings['bayes_dist_type'] - self.Bayesian_dict['pars_lower_bnds'] = input_dict['bounds']['lower'] - self.Bayesian_dict['pars_upper_bnds'] = input_dict['bounds']['upper'] + self.Bayesian_dict['rate_constants_initial_guess'] = self.x0; print("line 290", self.Bayesian_dict['pars_initial_guess']) #we **only** want the ones being changed. I assume you have this, but if not, there is a helper function I made called get_varying_rate_vals_and_bnds in CheKiPEUQ_from_Frhodo.py + self.Bayesian_dict['rate_constants_lower_bnds'] = input_dict['bounds']['lower']; print ("line 293", input_dict['bounds']['lower']) #we **only** want the ones being changed. + self.Bayesian_dict['rate_constants_upper_bnds'] = input_dict['bounds']['upper']; print ("line 294", input_dict['bounds']['lower']) #we **only** want the ones being changed. + + self.Bayesian_dict['rate_constants_parameters_initial_guess'] = [] #To be filled #we **only** want the ones being changed. + self.Bayesian_dict['rate_constants_parameters_lower_bnds'] = [] #To be filled #we **only** want the ones being changed. + self.Bayesian_dict['rate_constants_parameters_upper_bnds'] = #To be filled #we **only** want the ones being changed. def __call__(self, s, optimizing=True): def append_output(output_dict, calc_resid_output): @@ -358,13 +366,30 @@ def append_output(output_dict, calc_resid_output): elif self.opt_settings['obj_fcn_type'] == 'Bayesian': # TODO: we should bring in x_bnds (the coefficent bounds) so that we can use the elementary step coefficients for Bayesian rather than the rate_val values. Bayesian_dict = self.Bayesian_dict - Bayesian_dict['simulation_function'] = np.concatenate(output_dict['obs_sim_interp'], axis=0) + ''' + Start of variables Travis needs to fill. + ''' + Bayesian_dict['rate_constants_current_guess'] = [] #we **only** want the ones that are being changed, not the full set. + Bayesian_dict['rate_constants_parameters_current_guess'] = [] #we **only** want the ones that are being changed, not the full set. + ''' + End of variables Travis needs to fill. + ''' + Bayesian_dict['last_obs_sim_interp'] = np.concatenate(output_dict['obs_sim_interp'], axis=0) Bayesian_dict['observed_data'] = np.concatenate(output_dict['obs_exp'], axis=0) Bayesian_dict['observed_data_lower_bounds'] = [] Bayesian_dict['observed_data_upper_bounds'] = [] + def get_last_obs_sim_interp(varying_rate_vals): #A. Savara added this. It needs to be here. + try: + last_obs_sim_interp = Bayesian_dict['last_obs_sim_interp'] + last_obs_sim_interp = np.array(shock['last_obs_sim_interp']).T + except: + print("this isline 207! There may be an error occurring!") + last_obs_sim_interp = None + return last_obs_sim_interp + Bayesian_dict['simulation_function'] = get_last_obs_sim_interp #using wrapper that just returns the last_obs_sim_interp if np.size(loss_resid) == 1: # optimizing single experiment - Bayesian_dict['weights_data'] = aggregate_weights + Bayesian_dict['weights_data'] = np.array(output_dict['aggregate_weights'], dtype=object) else: aggregate_weights = np.array(output_dict['aggregate_weights'], dtype=object) exp_loss_weights = loss_exp/generalized_loss_fcn(loss_resid) # comparison is between selected loss fcn and SSE (L2 loss) @@ -375,19 +400,38 @@ def append_output(output_dict, calc_resid_output): #for val in Bayesian_dict['weights_data']: # print(val) - ''' - - - - - ASHI THIS IS WHERE YOU'D DO YOUR VOODOO MAGIC AND ALTER obj_fcn - - - - + + #Start of getting Bayesian objecive_function_value + import optimize.CheKiPEUQ_from_Frhodo + print("line 392", x) + #Step 3 of Bayesian: create a CheKiPEUQ_PE_Object (this is a class object) + #NOTE: normally, the Bayesian object would be created earlier. However, we are using a non-standard application + #where the observed_data uncertainties might change with each simulation. + #To gain efficiency, we could cheat and se the feature get_responses_simulation_uncertainties function of CheKiPEUQ, + #or could create a new get_responses_observed_uncertainties function in CheKiPEUQ + #for now we will just create a new PE_object each time. + CheKiPEUQ_PE_object = optimize.CheKiPEUQ_from_Frhodo.load_into_CheKiPUEQ( + simulation_function= Bayesian_dict['simulation_function'], + observed_data= Bayesian_dict['observed_data'], + pars_initial_guess = Bayesian_dict['pars_initial_guess'], #this is assigned in the "__init__" function up above. + pars_lower_bnds = Bayesian_dict['pars_lower_bnds'], #this is assigned in the "__init__" function up above. + pars_upper_bnds = Bayesian_dict['pars_upper_bnds'], #this is assigned in the "__init__" function up above. + observed_data_lower_bounds= Bayesian_dict['observed_data_lower_bounds'], + observed_data_upper_bounds= Bayesian_dict['observed_data_upper_bounds'], + weights_data= Bayesian_dict['weights_data'], + pars_uncertainty_distribution= Bayesian_dict['pars_uncertainty_distribution']) #this is assigned in the "__init__" function up above. + #Step 4 of Bayesian: call a function to get the posterior density which will be used as the objective function. + #We need to provide the current values of the varying_rate_vals to feed into the function. + print("line 406", varying_rate_vals_indices) + + + varying_rate_vals = np.array(output_dict['shock']['rate_val'])[list(varying_rate_vals_indices)] #when extracting a list of multiple indices, instead of array[index] one use array[[indices]] + log_posterior_density = optimize.CheKiPEUQ_from_Frhodo.get_log_posterior_density(CheKiPEUQ_PE_object, varying_rate_vals) + #Step 5 of Bayesian: return the objective function and any other metrics desired. + obj_fcn = -1*log_posterior_density #need neg_logP because minimizing. + #End of getting Bayesian objecive_function_value - ''' # For updating self.i += 1