popt.update_schemes.genopt
1# External imports 2import numpy as np 3from numpy import linalg as la 4import time 5 6# Internal imports 7from popt.misc_tools import optim_tools as ot 8from popt.loop.optimize import Optimize 9import popt.update_schemes.optimizers as opt 10from popt.update_schemes.cma import CMA 11 12 13class GenOpt(Optimize): 14 15 def __init__(self, fun, x, args, jac, jac_mut, corr_adapt=None, bounds=None, **options): 16 17 """ 18 Parameters 19 ---------- 20 fun: callable 21 objective function 22 23 x: ndarray 24 Initial state 25 26 args: tuple 27 Initial covariance 28 29 jac: callable 30 Gradient function 31 32 jac_mut: callable 33 Mutation gradient function 34 35 corr_adapt : callable 36 Function for correalation matrix adaption 37 38 bounds: list, optional 39 (min, max) pairs for each element in x. None is used to specify no bound. 40 41 options: dict 42 Optimization options 43 """ 44 45 # init PETEnsemble 46 super(GenOpt, self).__init__(**options) 47 48 def __set__variable(var_name=None, defalut=None): 49 if var_name in options: 50 return options[var_name] 51 else: 52 return defalut 53 54 # Set input as class variables 55 self.options = options # options 56 self.fun = fun # objective function 57 self.jac = jac # gradient function 58 self.jac_mut = jac_mut # mutation function 59 self.corr_adapt = corr_adapt # correlation adaption function 60 self.bounds = bounds # parameter bounds 61 self.mean_state = x # initial mean state 62 self.theta = args[0] # initial theta and correlation 63 self.corr = args[1] # inital correlation 64 65 # Set other optimization parameters 66 self.obj_func_tol = __set__variable('obj_func_tol', 1e-6) 67 self.alpha = __set__variable('alpha', 0.1) 68 self.alpha_theta = __set__variable('alpha_theta', 0.1) 69 self.alpha_corr = __set__variable('alpha_theta', 0.1) 70 self.beta = __set__variable('beta', 0.0) # this is stored in the optimizer class 71 self.nesterov = __set__variable('nesterov', False) # use Nesterov acceleration if value is true 72 self.alpha_iter_max = __set__variable('alpha_maxiter', 5) 73 self.max_resample = __set__variable('resample', 0) 74 self.normalize = __set__variable('normalize', True) 75 self.cov_factor = __set__variable('cov_factor', 0.5) 76 77 # Initialize other variables 78 self.state_step = 0 # state step 79 self.theta_step = 0 # covariance step 80 81 # Calculate objective function of startpoint 82 if not self.restart: 83 self.start_time = time.perf_counter() 84 self.obj_func_values = self.fun(self.mean_state) 85 self.nfev += 1 86 self.optimize_result = ot.get_optimize_result(self) 87 ot.save_optimize_results(self.optimize_result) 88 if self.logger is not None: 89 self.logger.info(' Running optimization...') 90 info_str = ' {:<10} {:<10} {:<15} {:<15} {:<10} {:<10} {:<10} {:<10} '.format('iter', 91 'alpha_iter', 92 'obj_func', 93 'step-size', 94 'alpha0', 95 'beta0', 96 'max corr', 97 'min_corr') 98 self.logger.info(info_str) 99 self.logger.info(' {:<21} {:<15.4e}'.format(self.iteration, 100 round(np.mean(self.obj_func_values),4))) 101 102 # Initialize optimizer 103 optimizer = __set__variable('optimizer', 'GA') 104 if optimizer == 'GA': 105 self.optimizer = opt.GradientAscent(self.alpha, self.beta) 106 elif optimizer == 'Adam': 107 self.optimizer = opt.Adam(self.alpha, self.beta) 108 109 # The GenOpt class self-ignites, and it is possible to send the EnOpt class as a callale method to scipy.minimize 110 self.run_loop() # run_loop resides in the Optimization class (super) 111 112 def calc_update(self): 113 """ 114 Update using steepest descent method with ensemble gradients 115 """ 116 117 # Initialize variables for this step 118 improvement = False 119 success = False 120 resampling_iter = 0 121 122 while improvement is False: # resampling loop 123 124 # Shrink covariance each time we try resampling 125 shrink = self.cov_factor ** resampling_iter 126 127 # Calculate gradient 128 if self.nesterov: 129 gradient = self.jac(self.mean_state + self.beta*self.state_step, 130 self.theta + self.beta*self.theta_step, self.corr) 131 else: 132 gradient = self.jac(self.mean_state, self.theta, self.corr) 133 self.njev += 1 134 135 # Compute the mutation gradient 136 gradient_theta, en_matrices = self.jac_mut(return_ensembles=True) 137 if self.normalize: 138 gradient /= np.maximum(la.norm(gradient, np.inf), 1e-12) # scale the gradient with inf-norm 139 gradient_theta /= np.maximum(la.norm(gradient_theta, np.inf), 1e-12) # scale the mutation with inf-norm 140 141 # Initialize for this step 142 alpha_iter = 0 143 144 while improvement is False: # backtracking loop 145 146 new_state, new_step = self.optimizer.apply_update(self.mean_state, gradient, iter=self.iteration) 147 new_state = ot.clip_state(new_state, self.bounds) 148 149 # Calculate new objective function 150 new_func_values = self.fun(new_state) 151 self.nfev += 1 152 153 if np.mean(self.obj_func_values) - np.mean(new_func_values) > self.obj_func_tol: 154 155 # Update objective function values and state 156 self.obj_func_values = new_func_values 157 self.mean_state = new_state 158 self.state_step = new_step 159 self.alpha = self.optimizer.get_step_size() 160 161 # Update theta (currently we don't apply backtracking for theta) 162 self.theta_step = self.beta*self.theta_step - self.alpha_theta*gradient_theta 163 self.theta = self.theta + self.theta_step 164 165 # update correlation matrix 166 if isinstance(self.corr_adapt, CMA): 167 enZ = en_matrices['gaussian'] 168 enJ = en_matrices['objective'] 169 self.corr = self.corr_adapt(cov = self.corr, 170 step = new_step/self.alpha, 171 X = enZ.T, 172 J = np.squeeze(enJ)) 173 174 elif callable(self.corr_adapt): 175 self.corr = self.corr - self.alpha_corr*self.corr_adapt() 176 177 # Write logging info 178 if self.logger is not None: 179 corr_max = round(np.max(self.corr-np.eye(self.corr.shape[0])), 3) 180 corr_min = round(np.min(self.corr), 3) 181 info_str_iter = ' {:<10} {:<10} {:<15.4e} {:<15} {:<10} {:<10} {:<10} {:<10}'.\ 182 format(self.iteration, 183 alpha_iter, 184 round(np.mean(self.obj_func_values),4), 185 self.alpha, 186 round(self.theta[0, 0],2), 187 round(self.theta[0, 1],2), 188 corr_max, 189 corr_min) 190 191 self.logger.info(info_str_iter) 192 193 # Update step size in the one-dimensional case 194 if new_state.size == 1: 195 self.optimizer.step_size /= 2 196 197 # Iteration was a success 198 improvement = True 199 success = True 200 self.optimizer.restore_parameters() 201 202 # Save variables defined in savedata keyword. 203 self.optimize_result = ot.get_optimize_result(self) 204 ot.save_optimize_results(self.optimize_result) 205 206 # Update iteration counter if iteration was successful and save current state 207 self.iteration += 1 208 209 else: 210 211 # If we do not have a reduction in the objective function, we reduce the step limiter 212 if alpha_iter < self.alpha_iter_max: 213 self.optimizer.apply_backtracking() # decrease alpha 214 alpha_iter += 1 215 elif (resampling_iter < self.max_resample and 216 np.mean(new_func_values) - np.mean(self.obj_func_values) > 0): # update gradient 217 resampling_iter += 1 218 self.optimizer.restore_parameters() 219 break 220 else: 221 success = False 222 return success 223 224 return success
14class GenOpt(Optimize): 15 16 def __init__(self, fun, x, args, jac, jac_mut, corr_adapt=None, bounds=None, **options): 17 18 """ 19 Parameters 20 ---------- 21 fun: callable 22 objective function 23 24 x: ndarray 25 Initial state 26 27 args: tuple 28 Initial covariance 29 30 jac: callable 31 Gradient function 32 33 jac_mut: callable 34 Mutation gradient function 35 36 corr_adapt : callable 37 Function for correalation matrix adaption 38 39 bounds: list, optional 40 (min, max) pairs for each element in x. None is used to specify no bound. 41 42 options: dict 43 Optimization options 44 """ 45 46 # init PETEnsemble 47 super(GenOpt, self).__init__(**options) 48 49 def __set__variable(var_name=None, defalut=None): 50 if var_name in options: 51 return options[var_name] 52 else: 53 return defalut 54 55 # Set input as class variables 56 self.options = options # options 57 self.fun = fun # objective function 58 self.jac = jac # gradient function 59 self.jac_mut = jac_mut # mutation function 60 self.corr_adapt = corr_adapt # correlation adaption function 61 self.bounds = bounds # parameter bounds 62 self.mean_state = x # initial mean state 63 self.theta = args[0] # initial theta and correlation 64 self.corr = args[1] # inital correlation 65 66 # Set other optimization parameters 67 self.obj_func_tol = __set__variable('obj_func_tol', 1e-6) 68 self.alpha = __set__variable('alpha', 0.1) 69 self.alpha_theta = __set__variable('alpha_theta', 0.1) 70 self.alpha_corr = __set__variable('alpha_theta', 0.1) 71 self.beta = __set__variable('beta', 0.0) # this is stored in the optimizer class 72 self.nesterov = __set__variable('nesterov', False) # use Nesterov acceleration if value is true 73 self.alpha_iter_max = __set__variable('alpha_maxiter', 5) 74 self.max_resample = __set__variable('resample', 0) 75 self.normalize = __set__variable('normalize', True) 76 self.cov_factor = __set__variable('cov_factor', 0.5) 77 78 # Initialize other variables 79 self.state_step = 0 # state step 80 self.theta_step = 0 # covariance step 81 82 # Calculate objective function of startpoint 83 if not self.restart: 84 self.start_time = time.perf_counter() 85 self.obj_func_values = self.fun(self.mean_state) 86 self.nfev += 1 87 self.optimize_result = ot.get_optimize_result(self) 88 ot.save_optimize_results(self.optimize_result) 89 if self.logger is not None: 90 self.logger.info(' Running optimization...') 91 info_str = ' {:<10} {:<10} {:<15} {:<15} {:<10} {:<10} {:<10} {:<10} '.format('iter', 92 'alpha_iter', 93 'obj_func', 94 'step-size', 95 'alpha0', 96 'beta0', 97 'max corr', 98 'min_corr') 99 self.logger.info(info_str) 100 self.logger.info(' {:<21} {:<15.4e}'.format(self.iteration, 101 round(np.mean(self.obj_func_values),4))) 102 103 # Initialize optimizer 104 optimizer = __set__variable('optimizer', 'GA') 105 if optimizer == 'GA': 106 self.optimizer = opt.GradientAscent(self.alpha, self.beta) 107 elif optimizer == 'Adam': 108 self.optimizer = opt.Adam(self.alpha, self.beta) 109 110 # The GenOpt class self-ignites, and it is possible to send the EnOpt class as a callale method to scipy.minimize 111 self.run_loop() # run_loop resides in the Optimization class (super) 112 113 def calc_update(self): 114 """ 115 Update using steepest descent method with ensemble gradients 116 """ 117 118 # Initialize variables for this step 119 improvement = False 120 success = False 121 resampling_iter = 0 122 123 while improvement is False: # resampling loop 124 125 # Shrink covariance each time we try resampling 126 shrink = self.cov_factor ** resampling_iter 127 128 # Calculate gradient 129 if self.nesterov: 130 gradient = self.jac(self.mean_state + self.beta*self.state_step, 131 self.theta + self.beta*self.theta_step, self.corr) 132 else: 133 gradient = self.jac(self.mean_state, self.theta, self.corr) 134 self.njev += 1 135 136 # Compute the mutation gradient 137 gradient_theta, en_matrices = self.jac_mut(return_ensembles=True) 138 if self.normalize: 139 gradient /= np.maximum(la.norm(gradient, np.inf), 1e-12) # scale the gradient with inf-norm 140 gradient_theta /= np.maximum(la.norm(gradient_theta, np.inf), 1e-12) # scale the mutation with inf-norm 141 142 # Initialize for this step 143 alpha_iter = 0 144 145 while improvement is False: # backtracking loop 146 147 new_state, new_step = self.optimizer.apply_update(self.mean_state, gradient, iter=self.iteration) 148 new_state = ot.clip_state(new_state, self.bounds) 149 150 # Calculate new objective function 151 new_func_values = self.fun(new_state) 152 self.nfev += 1 153 154 if np.mean(self.obj_func_values) - np.mean(new_func_values) > self.obj_func_tol: 155 156 # Update objective function values and state 157 self.obj_func_values = new_func_values 158 self.mean_state = new_state 159 self.state_step = new_step 160 self.alpha = self.optimizer.get_step_size() 161 162 # Update theta (currently we don't apply backtracking for theta) 163 self.theta_step = self.beta*self.theta_step - self.alpha_theta*gradient_theta 164 self.theta = self.theta + self.theta_step 165 166 # update correlation matrix 167 if isinstance(self.corr_adapt, CMA): 168 enZ = en_matrices['gaussian'] 169 enJ = en_matrices['objective'] 170 self.corr = self.corr_adapt(cov = self.corr, 171 step = new_step/self.alpha, 172 X = enZ.T, 173 J = np.squeeze(enJ)) 174 175 elif callable(self.corr_adapt): 176 self.corr = self.corr - self.alpha_corr*self.corr_adapt() 177 178 # Write logging info 179 if self.logger is not None: 180 corr_max = round(np.max(self.corr-np.eye(self.corr.shape[0])), 3) 181 corr_min = round(np.min(self.corr), 3) 182 info_str_iter = ' {:<10} {:<10} {:<15.4e} {:<15} {:<10} {:<10} {:<10} {:<10}'.\ 183 format(self.iteration, 184 alpha_iter, 185 round(np.mean(self.obj_func_values),4), 186 self.alpha, 187 round(self.theta[0, 0],2), 188 round(self.theta[0, 1],2), 189 corr_max, 190 corr_min) 191 192 self.logger.info(info_str_iter) 193 194 # Update step size in the one-dimensional case 195 if new_state.size == 1: 196 self.optimizer.step_size /= 2 197 198 # Iteration was a success 199 improvement = True 200 success = True 201 self.optimizer.restore_parameters() 202 203 # Save variables defined in savedata keyword. 204 self.optimize_result = ot.get_optimize_result(self) 205 ot.save_optimize_results(self.optimize_result) 206 207 # Update iteration counter if iteration was successful and save current state 208 self.iteration += 1 209 210 else: 211 212 # If we do not have a reduction in the objective function, we reduce the step limiter 213 if alpha_iter < self.alpha_iter_max: 214 self.optimizer.apply_backtracking() # decrease alpha 215 alpha_iter += 1 216 elif (resampling_iter < self.max_resample and 217 np.mean(new_func_values) - np.mean(self.obj_func_values) > 0): # update gradient 218 resampling_iter += 1 219 self.optimizer.restore_parameters() 220 break 221 else: 222 success = False 223 return success 224 225 return success
Class for ensemble optimization algorithms. These are classified by calculating the sensitivity or gradient using ensemble instead of classical derivatives. The loop is else as a classic optimization loop: a state (or control variable) will be iterated upon using an algorithm defined in the update_scheme package.
Attributes
- logger (Logger): Print output to screen and log-file
- pickle_restart_file (str): Save name for pickle dump/load
- optimize_result (OptimizeResult): Dictionary with results for the current iteration
- iteration (int): Iteration index
- max_iter (int): Max number of iterations
- restart (bool): Restart flag
- restartsave (bool): Save restart information flag
Methods
run_loop() The main optimization loop
save() Save restart file
load() Load restart file
calc_update() Empty dummy function, actual functionality must be defined by the subclasses
GenOpt(fun, x, args, jac, jac_mut, corr_adapt=None, bounds=None, **options)
16 def __init__(self, fun, x, args, jac, jac_mut, corr_adapt=None, bounds=None, **options): 17 18 """ 19 Parameters 20 ---------- 21 fun: callable 22 objective function 23 24 x: ndarray 25 Initial state 26 27 args: tuple 28 Initial covariance 29 30 jac: callable 31 Gradient function 32 33 jac_mut: callable 34 Mutation gradient function 35 36 corr_adapt : callable 37 Function for correalation matrix adaption 38 39 bounds: list, optional 40 (min, max) pairs for each element in x. None is used to specify no bound. 41 42 options: dict 43 Optimization options 44 """ 45 46 # init PETEnsemble 47 super(GenOpt, self).__init__(**options) 48 49 def __set__variable(var_name=None, defalut=None): 50 if var_name in options: 51 return options[var_name] 52 else: 53 return defalut 54 55 # Set input as class variables 56 self.options = options # options 57 self.fun = fun # objective function 58 self.jac = jac # gradient function 59 self.jac_mut = jac_mut # mutation function 60 self.corr_adapt = corr_adapt # correlation adaption function 61 self.bounds = bounds # parameter bounds 62 self.mean_state = x # initial mean state 63 self.theta = args[0] # initial theta and correlation 64 self.corr = args[1] # inital correlation 65 66 # Set other optimization parameters 67 self.obj_func_tol = __set__variable('obj_func_tol', 1e-6) 68 self.alpha = __set__variable('alpha', 0.1) 69 self.alpha_theta = __set__variable('alpha_theta', 0.1) 70 self.alpha_corr = __set__variable('alpha_theta', 0.1) 71 self.beta = __set__variable('beta', 0.0) # this is stored in the optimizer class 72 self.nesterov = __set__variable('nesterov', False) # use Nesterov acceleration if value is true 73 self.alpha_iter_max = __set__variable('alpha_maxiter', 5) 74 self.max_resample = __set__variable('resample', 0) 75 self.normalize = __set__variable('normalize', True) 76 self.cov_factor = __set__variable('cov_factor', 0.5) 77 78 # Initialize other variables 79 self.state_step = 0 # state step 80 self.theta_step = 0 # covariance step 81 82 # Calculate objective function of startpoint 83 if not self.restart: 84 self.start_time = time.perf_counter() 85 self.obj_func_values = self.fun(self.mean_state) 86 self.nfev += 1 87 self.optimize_result = ot.get_optimize_result(self) 88 ot.save_optimize_results(self.optimize_result) 89 if self.logger is not None: 90 self.logger.info(' Running optimization...') 91 info_str = ' {:<10} {:<10} {:<15} {:<15} {:<10} {:<10} {:<10} {:<10} '.format('iter', 92 'alpha_iter', 93 'obj_func', 94 'step-size', 95 'alpha0', 96 'beta0', 97 'max corr', 98 'min_corr') 99 self.logger.info(info_str) 100 self.logger.info(' {:<21} {:<15.4e}'.format(self.iteration, 101 round(np.mean(self.obj_func_values),4))) 102 103 # Initialize optimizer 104 optimizer = __set__variable('optimizer', 'GA') 105 if optimizer == 'GA': 106 self.optimizer = opt.GradientAscent(self.alpha, self.beta) 107 elif optimizer == 'Adam': 108 self.optimizer = opt.Adam(self.alpha, self.beta) 109 110 # The GenOpt class self-ignites, and it is possible to send the EnOpt class as a callale method to scipy.minimize 111 self.run_loop() # run_loop resides in the Optimization class (super)
Parameters
- fun (callable): objective function
- x (ndarray): Initial state
- args (tuple): Initial covariance
- jac (callable): Gradient function
- jac_mut (callable): Mutation gradient function
- corr_adapt (callable): Function for correalation matrix adaption
- bounds (list, optional): (min, max) pairs for each element in x. None is used to specify no bound.
- options (dict): Optimization options
def
calc_update(self):
113 def calc_update(self): 114 """ 115 Update using steepest descent method with ensemble gradients 116 """ 117 118 # Initialize variables for this step 119 improvement = False 120 success = False 121 resampling_iter = 0 122 123 while improvement is False: # resampling loop 124 125 # Shrink covariance each time we try resampling 126 shrink = self.cov_factor ** resampling_iter 127 128 # Calculate gradient 129 if self.nesterov: 130 gradient = self.jac(self.mean_state + self.beta*self.state_step, 131 self.theta + self.beta*self.theta_step, self.corr) 132 else: 133 gradient = self.jac(self.mean_state, self.theta, self.corr) 134 self.njev += 1 135 136 # Compute the mutation gradient 137 gradient_theta, en_matrices = self.jac_mut(return_ensembles=True) 138 if self.normalize: 139 gradient /= np.maximum(la.norm(gradient, np.inf), 1e-12) # scale the gradient with inf-norm 140 gradient_theta /= np.maximum(la.norm(gradient_theta, np.inf), 1e-12) # scale the mutation with inf-norm 141 142 # Initialize for this step 143 alpha_iter = 0 144 145 while improvement is False: # backtracking loop 146 147 new_state, new_step = self.optimizer.apply_update(self.mean_state, gradient, iter=self.iteration) 148 new_state = ot.clip_state(new_state, self.bounds) 149 150 # Calculate new objective function 151 new_func_values = self.fun(new_state) 152 self.nfev += 1 153 154 if np.mean(self.obj_func_values) - np.mean(new_func_values) > self.obj_func_tol: 155 156 # Update objective function values and state 157 self.obj_func_values = new_func_values 158 self.mean_state = new_state 159 self.state_step = new_step 160 self.alpha = self.optimizer.get_step_size() 161 162 # Update theta (currently we don't apply backtracking for theta) 163 self.theta_step = self.beta*self.theta_step - self.alpha_theta*gradient_theta 164 self.theta = self.theta + self.theta_step 165 166 # update correlation matrix 167 if isinstance(self.corr_adapt, CMA): 168 enZ = en_matrices['gaussian'] 169 enJ = en_matrices['objective'] 170 self.corr = self.corr_adapt(cov = self.corr, 171 step = new_step/self.alpha, 172 X = enZ.T, 173 J = np.squeeze(enJ)) 174 175 elif callable(self.corr_adapt): 176 self.corr = self.corr - self.alpha_corr*self.corr_adapt() 177 178 # Write logging info 179 if self.logger is not None: 180 corr_max = round(np.max(self.corr-np.eye(self.corr.shape[0])), 3) 181 corr_min = round(np.min(self.corr), 3) 182 info_str_iter = ' {:<10} {:<10} {:<15.4e} {:<15} {:<10} {:<10} {:<10} {:<10}'.\ 183 format(self.iteration, 184 alpha_iter, 185 round(np.mean(self.obj_func_values),4), 186 self.alpha, 187 round(self.theta[0, 0],2), 188 round(self.theta[0, 1],2), 189 corr_max, 190 corr_min) 191 192 self.logger.info(info_str_iter) 193 194 # Update step size in the one-dimensional case 195 if new_state.size == 1: 196 self.optimizer.step_size /= 2 197 198 # Iteration was a success 199 improvement = True 200 success = True 201 self.optimizer.restore_parameters() 202 203 # Save variables defined in savedata keyword. 204 self.optimize_result = ot.get_optimize_result(self) 205 ot.save_optimize_results(self.optimize_result) 206 207 # Update iteration counter if iteration was successful and save current state 208 self.iteration += 1 209 210 else: 211 212 # If we do not have a reduction in the objective function, we reduce the step limiter 213 if alpha_iter < self.alpha_iter_max: 214 self.optimizer.apply_backtracking() # decrease alpha 215 alpha_iter += 1 216 elif (resampling_iter < self.max_resample and 217 np.mean(new_func_values) - np.mean(self.obj_func_values) > 0): # update gradient 218 resampling_iter += 1 219 self.optimizer.restore_parameters() 220 break 221 else: 222 success = False 223 return success 224 225 return success
Update using steepest descent method with ensemble gradients