"""ParticleSwarm optimizer (the ``pso`` fit type).
Extracted from the algorithms package (M1 Step 4). Subclasses Algorithm and
inherits its run loop + execution seam; it makes no core.* call of its own.
"""
from ..base import Algorithm
from ...config_schema import PyBNFConfigModel
from ...pset import PSet
from ...printing import print1, print2
from ...registry import register_fit_type
import logging
import numpy as np
import re
# Preserve the original module logger name so log records keep the
# 'pybnf.algorithms' channel.
logger = logging.getLogger('pybnf.algorithms')
[docs]
class PSOConfig(PyBNFConfigModel):
"""Particle-swarm config fields, co-located with the method (ADR-0002,
ADR-0006). These defaults previously lived in ``config_schema.GlobalConfig``;
moving them here makes ``pso``'s config knowledge travel with its algorithm
class. The values are byte-identical to the old global defaults -- the
int-literal ``adaptive_n_max = 30`` stays an int (pydantic does not validate
defaults), matching the golden-config net.
``particle_weight_final`` is intentionally NOT here: it is not a default but a
runtime fallback (set to ``particle_weight`` when absent, in ``__init__``
below), so it stays a pass-through extra.
"""
particle_weight: float = 0.7
adaptive_n_max: float = 30
adaptive_n_stop: float = float('inf')
adaptive_abs_tol: float = 0.0
adaptive_rel_tol: float = 0.0
cognitive: float = 1.5
social: float = 1.5
v_stop: float = 0.0
# particle_weight_final is a runtime fallback (-> particle_weight when absent, set
# in __init__), not a schema default, but it IS a valid pso key (#401).
RUNTIME_KEYS = frozenset({'particle_weight_final'})
[docs]
@register_fit_type('pso', family='optimizer', display_name='Particle Swarm Optimization',
schema=PSOConfig)
class ParticleSwarm(Algorithm):
"""
Implements particle swarm optimization.
The implementation roughly follows Moraes et al 2015, although is reorganized to better suit PyBNF's format.
Note the global convergence criterion discussed in that paper is not used (would require too long a
computation), and instead uses ????
"""
def __init__(self, config):
# Former params that are now part of the config
# variable_list, num_particles, max_evals, cognitive=1.5, social=1.5, w0=1.,
# wf=0.1, nmax=30, n_stop=np.inf, absolute_tol=0., relative_tol=0.)
"""
Initial configuration of particle swarm optimizer
:param conf_dict: The fitting configuration
:type conf_dict: Configuration
The config should contain the following definitions:
population_size - Number of particles in the swarm
max_iterations - Maximum number of iterations. More precisely, the max number of simulations run is this times
the population size.
cognitive - Acceleration toward the particle's own best
social - Acceleration toward the global best
particle_weight - Inertia weight of the particle (default 1)
The following config parameters relate to the complicated method presented is Moraes et al for adjusting the
inertia weight as you go. These are optional, and this feature will be disabled (by setting
particle_weight_final = particle_weight) if these are not included.
It remains to be seen whether this method is at all useful for our applications.
particle_weight_final - Inertia weight at the end of the simulation
adaptive_n_max - Controls how quickly we approach wf - After nmax "unproductive" iterations, we are halfway from
w0 to wf
adaptive_n_stop - nd the entire run if we have had this many "unproductive" iterations (should be more than
adaptive_n_max)
adaptive_abs_tol - Tolerance for determining if an iteration was "unproductive". A run is unproductive if the
change in global_best is less than absolute_tol + relative_tol * global_best
adaptive_rel_tol - Tolerance 2 for determining if an iteration was "unproductive" (see above)
"""
super().__init__(config)
# This default value gets special treatment because if missing, it should take the value of particle_weight,
# disabling the adaptive weight change entirely.
if 'particle_weight_final' not in self.config.config:
self.config.config['particle_weight_final'] = self.config.config['particle_weight']
# Save config parameters
self.c1 = self.config.config['cognitive']
self.c2 = self.config.config['social']
self.max_evals = self.config.config['population_size'] * self.config.config['max_iterations']
self.output_every = self.config.config['population_size'] * self.config.config['output_every']
self.num_particles = self.config.config['population_size']
# Todo: Nice error message if a required key is missing
self.w0 = self.config.config['particle_weight']
self.wf = self.config.config['particle_weight_final']
self.nmax = self.config.config['adaptive_n_max']
self.n_stop = self.config.config['adaptive_n_stop']
self.absolute_tol = self.config.config['adaptive_abs_tol']
self.relative_tol = self.config.config['adaptive_rel_tol']
self.nv = 0 # Counter that controls the current weight. Counts number of "unproductive" iterations.
self.num_evals = 0 # Counter for the total number of results received
# Initialize storage for the swarm data
self.swarm = [] # List of lists of the form [PSet, velocity]. Velocity is stored as a dict with the same keys
# as PSet
self.pset_map = dict() # Maps each PSet to it s particle number, for easy lookup.
# One independent [PSet, objective] slot per particle. A list-multiply
# ([[None, inf]] * n) would alias one inner list across every particle --
# harmless today (got_result reassigns the whole slot) but a foot-gun if
# anyone ever mutates a slot in place.
self.bests = [[None, np.inf] for _ in range(self.num_particles)] # The best result for each particle
self.global_best = [None, np.inf] # The best result for the whole swarm
self.last_best = np.inf
[docs]
def reset(self, bootstrap=None):
super().reset(bootstrap)
self.nv = 0
self.num_evals = 0
self.swarm = []
self.pset_map = dict()
self.bests = [[None, np.inf] for _ in range(self.num_particles)]
self.global_best = [None, np.inf]
self.last_best = np.inf
[docs]
def start_run(self):
"""
Start the run by initializing n particles at random positions and velocities
:return:
"""
print2('Running Particle Swarm Optimization with %i particles for %i total simulations' %
(self.num_particles, self.max_evals))
if self.config.config['initialization'] == 'lh':
new_params_list = self.random_latin_hypercube_psets(self.num_particles)
else:
new_params_list = [self.random_pset() for i in range(self.num_particles)]
# ADR-0043 Phase 2: seed exactly one particle at the initial_value point (a no-op
# unless a parameter: record declares one); the rest stay random for diversity.
new_params_list[0] = self._seed_initial_value_pset(new_params_list[0])
for i in range(len(new_params_list)):
p = new_params_list[i]
p.name = 'iter0p%i' % i
# As suggested by Engelbrecht 2012, set all initial velocities to 0
new_velocity = dict({v.name: 0. for v in self.variables})
self.swarm.append([p, new_velocity])
self.pset_map[p] = len(self.swarm)-1 # Index of the newly added PSet.
return [particle[0] for particle in self.swarm]
[docs]
def got_result(self, res):
"""
Updates particle velocity and position after a simulation completes.
:param res: Result object containing the run PSet and the resulting Data.
:return:
"""
paramset = res.pset
score = res.score
self.num_evals += 1
if self.num_evals % self.num_particles == 0:
if (self.num_evals / self.num_particles) % 10 == 0:
print1('Completed %i of %i simulations' % (self.num_evals, self.max_evals))
else:
print2('Completed %i of %i simulations' % (self.num_evals, self.max_evals))
print2(f'Current best score: {self.global_best[1]:f}')
# End of one "pseudoflight", check if it was productive.
if (self.last_best != np.inf and
np.abs(self.last_best - self.global_best[1]) <
self.absolute_tol + self.relative_tol * self.last_best):
self.nv += 1
self.last_best = self.global_best[1]
# Check stop criterion
if self.config.config['v_stop'] > 0:
max_speed = max([abs(v) for p in self.swarm for v in p[1].values()])
if max_speed < self.config.config['v_stop']:
logger.info(f'Stopping particle swarm because the max speed is {max_speed}')
return 'STOP'
if self.num_evals % self.output_every == 0:
self.output_results()
p = self.pset_map.pop(paramset) # Particle number
# Update best scores if needed.
if score <= self.bests[p][1]:
self.bests[p] = [paramset, score]
if score <= self.global_best[1]:
self.global_best = [paramset, score]
# Update own position and velocity
# The order matters - updating velocity first seems to make the best use of our current info.
w = self.w0 + (self.wf - self.w0) * self.nv / (self.nv + self.nmax)
self.swarm[p][1] = \
{v.name:
w * self.swarm[p][1][v.name] +
self.c1 * self.rng.random() * self.bests[p][0].get_param(v.name).diff(self.swarm[p][0].get_param(v.name)) +
self.c2 * self.rng.random() * self.global_best[0].get_param(v.name).diff(self.swarm[p][0].get_param(v.name))
for v in self.variables}
# Manually check to determine if reflection occurred (i.e. attempted assigning of variable outside its bounds)
# If so, update based on reflection protocol and set velocity to 0
new_vars = []
for v in self.swarm[p][0]:
new_vars.append(v.add(self.swarm[p][1][v.name]))
# The would-be value after the velocity step, in sampling space u then
# mapped back via the parameter's scale -- to detect (without applying)
# an out-of-box move and zero the velocity (#412).
new_val = v.from_sampling_space(v.to_sampling_space(v.value) + self.swarm[p][1][v.name])
if new_val < v.lower_bound or v.upper_bound < new_val:
self.swarm[p][1][v.name] = 0.0
new_pset = PSet(new_vars)
self.swarm[p][0] = new_pset
# This will cause a crash if new_pset happens to be the same as an already running pset in pset_map.
# This could come up in practice if all parameters have hit a box constraint.
# As a simple workaround, perturb the parameters slightly
while new_pset in self.pset_map:
new_pset = PSet([v.add_rand(-1e-6, 1e-6, self.rng) for v in self.swarm[p][0]])
self.pset_map[new_pset] = p
# Set the new name: the old pset name is iter##p##
# Extract the iter number
iternum = int(re.search('iter([0-9]+)', paramset.name).groups()[0])
new_pset.name = 'iter%ip%i' % (iternum+1, p)
# Check for stopping criteria
if self.num_evals >= self.max_evals or self.nv >= self.n_stop:
return 'STOP'
return [new_pset]
[docs]
def add_iterations(self, n):
self.max_evals += n * self.config.config['population_size']