Processes 2019, 7, 163
the global optimization algorithm implemented in the bifurcation evolution tool, performs poorly
with large parameter value ranges, so selecting appropriate bounds is critical. To accommodate for
this, the automatic selection, which is the default setting in the tool, attempts to narrow the ranges
in a generalizable manner. First, the algorithm checks the current parameter value, p i , in the loaded
model, and, if the value is zero, the algorithm creates a range of parameter values from 1 × 10 −25 to 10,
approximating zero but prohibiting possible failure of the algorithm if the parameter appears in the
denominator of a rate law. For parameter values less than or equal to 10, the assigned parameter range
is from
p i
10 to 10. All parameter values greater than 10 receive a range from
p i
10 to 2p i . The automatic
assignment allows appropriate flexibility if the range of suitable parameter values for bifurcation is
unknown, but relies on the initial parameter value to provide a suitable estimate. If an appropriate
range is available for a given parameter through reliable experimental data, manual assignment may
be preferable, particularly if this assignment further narrows the range.
Within the differential evolution global optimization algorithm, parameter values are selected
from the assigned ranges using a random uniform distribution, such that it is equally likely to choose
any value within the assigned range. As this would greatly reduce the frequency of assigned parameter
values less than 1, and possibly prevent a parameter from occupying the ideal parameter space to
achieve the desired bifurcation behavior, the algorithm selects from a random uniform log distribution
of the parameter ranges for all parameters with an upper bound of 10. This ensures that values across
multiple orders of magnitude are equally likely to be selected.
2.3. Differential Evolution Algorithm
In the Pythonic approach to this tool, a simple implementation of the differential evolution
algorithm developed by Storn and Price was integrated within the bifurcation module, to perform a
search in parameter space for global minima of the objective function [13].
2.3.1. Initializing a Population
The algorithm begins by initializing a population of parameter vectors which represent solutions
to the objective function. Parameter space in this multi-dimensional optimization problem contains N
dimensions, where N is the number of parameters being optimized. These vectors are populated with
elements assigned randomly selected values within the predefined bounds specified for N i , a single
parameter, such that the members in the initial population occupy diverse regions of parameter
space. While Pythonic versions of this algorithm have been developed, the version available through
the scipy.optimize package, frequently used for similar optimization problems, does not allow unfit
members to be discarded from the starting population [14]. This makes optimization inefficient, slows
convergence and increases the likelihood that the algorithm will terminate before a sufficient minima is
reached. The current implementation of the algorithm discards all members with an objective function
evaluation above a predetermined threshold before evolving the population. This threshold value
coincides with penalty functions included within the bifurcation objective function so that parameter
vectors which do not reach steady state, or which have multiple eigenvalues approaching zero in both
the real and complex component, are discarded from the solution.
2.3.2. Recombination
During a single round of differential evolution, each member of the population undergoes
recombination to construct a trial vector. While iterating through each element, the trial vector is
populated with parameter values taken from the member at the current population index or from a
mutant vector. If a random number chosen from between zero and one is smaller than the crossover
probability, the trial vector receives the parameter value from the mutant vector, as long as the
parameter value remains within the acceptable range. Otherwise, the trial vector receives the element
from the current population member.
7
the global optimization algorithm implemented in the bifurcation evolution tool, performs poorly
with large parameter value ranges, so selecting appropriate bounds is critical. To accommodate for
this, the automatic selection, which is the default setting in the tool, attempts to narrow the ranges
in a generalizable manner. First, the algorithm checks the current parameter value, p i , in the loaded
model, and, if the value is zero, the algorithm creates a range of parameter values from 1 × 10 −25 to 10,
approximating zero but prohibiting possible failure of the algorithm if the parameter appears in the
denominator of a rate law. For parameter values less than or equal to 10, the assigned parameter range
is from
p i
10 to 10. All parameter values greater than 10 receive a range from
p i
10 to 2p i . The automatic
assignment allows appropriate flexibility if the range of suitable parameter values for bifurcation is
unknown, but relies on the initial parameter value to provide a suitable estimate. If an appropriate
range is available for a given parameter through reliable experimental data, manual assignment may
be preferable, particularly if this assignment further narrows the range.
Within the differential evolution global optimization algorithm, parameter values are selected
from the assigned ranges using a random uniform distribution, such that it is equally likely to choose
any value within the assigned range. As this would greatly reduce the frequency of assigned parameter
values less than 1, and possibly prevent a parameter from occupying the ideal parameter space to
achieve the desired bifurcation behavior, the algorithm selects from a random uniform log distribution
of the parameter ranges for all parameters with an upper bound of 10. This ensures that values across
multiple orders of magnitude are equally likely to be selected.
2.3. Differential Evolution Algorithm
In the Pythonic approach to this tool, a simple implementation of the differential evolution
algorithm developed by Storn and Price was integrated within the bifurcation module, to perform a
search in parameter space for global minima of the objective function [13].
2.3.1. Initializing a Population
The algorithm begins by initializing a population of parameter vectors which represent solutions
to the objective function. Parameter space in this multi-dimensional optimization problem contains N
dimensions, where N is the number of parameters being optimized. These vectors are populated with
elements assigned randomly selected values within the predefined bounds specified for N i , a single
parameter, such that the members in the initial population occupy diverse regions of parameter
space. While Pythonic versions of this algorithm have been developed, the version available through
the scipy.optimize package, frequently used for similar optimization problems, does not allow unfit
members to be discarded from the starting population [14]. This makes optimization inefficient, slows
convergence and increases the likelihood that the algorithm will terminate before a sufficient minima is
reached. The current implementation of the algorithm discards all members with an objective function
evaluation above a predetermined threshold before evolving the population. This threshold value
coincides with penalty functions included within the bifurcation objective function so that parameter
vectors which do not reach steady state, or which have multiple eigenvalues approaching zero in both
the real and complex component, are discarded from the solution.
2.3.2. Recombination
During a single round of differential evolution, each member of the population undergoes
recombination to construct a trial vector. While iterating through each element, the trial vector is
populated with parameter values taken from the member at the current population index or from a
mutant vector. If a random number chosen from between zero and one is smaller than the crossover
probability, the trial vector receives the parameter value from the mutant vector, as long as the
parameter value remains within the acceptable range. Otherwise, the trial vector receives the element
from the current population member.
7