# Surrogate-assisted distributed swarm optimisation for computationally
expensive models

Rohitash Chandra<EMAIL_ADDRESS>[ Yash Vardhan Sharma
<EMAIL_ADDRESS>Transitional Artificial Intelligence Research Group,
School of Mathematics and Statistics, University of New South Wales, Sydney,
NSW 2052, Australia  Mechanical and Industrial Engineering Department, Indian
Institute of Technology Roorkee, India

###### Abstract

Advances in parallel and distributed computing have enabled efficient
implementation of distributed swarm and evolutionary algorithms for complex
and computationally expensive models. Evolutionary algorithms provide
gradient-free optimisation which is beneficial for models that do not have
such information available, for instance geoscientific landscape evolution
models. However, such models are so computationally expensive that even
distributed swarm and evolutionary algorithms with the power of parallel
computing struggle. We need to incorporate efficient strategies such as
surrogate assisted optimisation that further improves their performance;
however this becomes a challenge given parallel processing and inter-process
communication for implementing surrogate training and prediction. In this
paper, we implement surrogate-based estimation of fitness evaluation in
distributed swarm optimisation over a parallel computing architecture. Our
results demonstrate very promising results for benchmark functions and
geoscientific landscape evolution models. We obtain reduction in
computationally time while retaining optimisation solution accuracy though the
use of surrogates in a parallel computing environment.

###### keywords:

Surrogate-assisted optimization, parallel computing, distributed swarm
algorithms, landscape evolution models, geoscientific models

††journal:

url]rohitash-chandra.github.io

##  1 Introduction

Evolutionary algorithms are loosely motivated by the theory of evolution where
species are represented by individuals in a population that compete and
collaborate with each other, producing offspring over generations that to
improve quality given by fitness measure [1, 2]. Particle swarm optimisation
(PSO) on the other hand is motivated by flocking behavior of birds or swarms
that is represented by a population of particles (individuals) that compete
and collaborate over time [3, 4, 5]. Evolutionary and swarm optimisation
methods have been prominent in a number of areas such as real-parameter global
optimization, combinatorial optimization and scheduling, and machine learning
[6, 7, 8, 3]. Research in evolutionary and swarm optimisation has focused on
different ways to create new solutions with mechanisms that are heuristic in
nature and hence different variants are available [9, 4]. A major challenge
has been in applying them in large scale or computationally expensive
optimisation problems that require thousands of function evaluations where a
single function evaluation can take minutes to hours or even days [10]. An
example of an expensive function is a geo-scientific model for landscape
evolution [11], and deep learning models for big data problems [12].
Computationally expensive optimisation can be addressed with distributed and
parallel computing with evolutionary and swarm optimisation [13, 14, 15];
however, we need efficient strategies for representing the problem.

Surrogate assisted optimisation [16, 17, 18, 19] provides a remedy for
expensive models with use of statistical and machine learning models to
provide a low computational replicate of the actual model. The replica or
surrogate is developed by training from available data generated during
optimisation that features a set of inputs (new solutions) and corresponding
output (fitness) given by actual model. The method is also called Bayesian
optimisation where the surrogate model which is known as acquisition function
is typically a Gaussian process model [20]; however, neural networks are also
used. Evolutionary and swarm optimization methods have been used in surrogate
assisted and Bayesian optimization, and been prominent in fields of engine and
aerospace design [21, 22, 23, 16], robotics [20], experimental design [24],
and machine learning [25, 26].

Evolutionary algorithms provide gradient-free optimisation which is beneficial
for models such as those that do not have gradient information, for instance
geo-scientific landscape evolution models [27]. Some instances of such models
are so expensive that even distributed evolutionary algorithms with the power
of parallel computing would struggle. Hence, we need to incorporate efficient
strategies such as surrogate assisted optimisation that further improves their
performance but this becomes a challenge given parallel processing and inter-
process communication for implementing surrogate estimation.

Landscape evolution models are geo-scientific models that has been used to
reconstruct the evolution of the Earth’s landscape over millions of years [28,
29, 30]. These models guide geologists and climate scientists in better
understanding Earths geologic and climate history that can further also help
in foreseeing the distant future of the planet [30, 31]. These models use data
from geological observations such as bore-hole data and estimated landscape
topography millions of years ago and require climate conditions and geological
parameters that are not easily known. Hence, it becomes an optimisation
problem to estimate these parameters which has been tackled mostly with
Bayesian inference via Markov Chain Monte Carlo (MCMC) sampling in previous
studies that used parallel computing [27], and surrogate assisted Bayesian
inference [32]. Motivated by these studies, we bring the problem of landscape
evolution model to the evolutionary and swarm optimisation community, rather
than viewing it is an inference problem, to see it as an optimisation problem.

In this paper, we implement surrogate-based optimisation framework via swarm
optimisation over a parallel computing architecture. We apply the framework
for benchmark optimisation functions and a selected landscape evolution model.
We investigate performance measures such as accuracy of surrogate prediction
given different types of problems, in terms of dimension and fitness
landscape.

The rest of the paper is organised as follows. Section 2 provides background
and related work, while Section 3 presents the proposed methodology. Section 4
presents experiments and results and Section 5 concludes the paper with
discussion for future research.

##  2 Related work

###  2.1 Distributed and parallel evolutionary and swarm optimisation

The area of swam optimisation is a growing field with innovative and
multidisciplinary optimisation strategy focused in improving quality of
solutions. Lin et. al [33] enhanced genetic learning particle swarm
optimization enhances PSO by breeding better solutions to guide the motion of
particles with diversity enhancement by ring topology. Liu et. al [34]
modified particle swarm optimization using adaptive strategy to avoid the
premature convergence, stochastic and mainstream learning strategies, an
adaptive position updating strategy and terminal replacement mechanism for
complex optimization problems. Roshanzamir et. al [35] presented hierarchical
multi group particle swarm optimization with different task allocations
inspired by holonic multi agent systems to improve their solutions.

In parallel computing, problem using separate parallel processors which may
feature inter-process communication typically within a single machine [36].
Distributed computing on the other hand, considers several machines to break
down a major task and therefore, distributed computing may also use parallel
computing [37]. Distributed evolutionary algorithms and distributed swarms
[38] typically refer to population distributed and dimension distributed
models. Distributed evolutionary and swarm optimisation focuses in enhancing
diversity of new populations and also to provide better computational
efficiency especially in high high dimensional problems [39, 38].

Distributed and parallel optimisation has been steadily making progress. One
of the earlier attempts was by Hereford [40] who presented distributed PSO
suitable for a swarm consisting of a large number of small, mobile robots
where each robot calculates its new position based on its present position and
present measurement eliminating need for central agent and had limited inter-
robot communication. Abadlia et al [41] used island model concepts with PSO
algorithm to improve its diversity and convergence, where the particles are
distrusted into separate islands with exchange through a migration process.
The results indicated that the island model provides improved efficiency of
PSO performance. Decampos et. al presented [42] parallel PSO strategies using
Pareto dominance and decomposition based on multiple swarms for many-objective
optimization problems and investigated the impact of synchronous and
asynchronous communication strategies for the decomposition-based approach.
The results indicated that parallelization had a positive effect on the
convergence and diversity of the optimization process. Erdeljan et al [43]
presented a data model partitioning for parallelization of analytical power
calculations using distributed PSO algorithms where results show that
distributed PSO algorithm achieves significantly better results than the basic
PSO algorithm.

###  2.2 Surrogate-assisted optimization

Surrogate-assisted optimisation has been initially motivated by
computationally expensive engineering and design problems. Ong _et al._ [18]
proposed parallel evolutionary optimisation with application to aerodynamic
wing design where surrogate models used radial basis functions (RBF) networks.
Zhou _et al._ [19] utilised global and local surrogate models for evolutionary
optimization and improved performance. Lim _et al._ [44] accounted for
uncertainty in estimation with a generalised method to unify diverse surrogate
models. The search for an appropriate surrogate model is a major challenge
given different forms of error or fitness landscape given by the actual model.
Giunta _et al._ [45] compared surrogate models that used quadratic polynomial
models with Gaussian process regression models (kriging) and found that the
former gave more accurate estimation for the given optimisation problems. Jin
_et al._ [46] compared several surrogate models where they reported radial
basis functions as one of the best for scalability and robustness and kriging
to be computationally expensive. There exists some comprehensive reviews about
different categories of surrogate-assisted optimisation for different types of
theoretical problems [17, 47] and applications in Earth science and water
resource modelling [48].

Swarm optimisation is becoming prominent in surrogate assisted optimisation.
Yu et. al [49] presented surrogate-assisted hierarchical particle swarm
optimization with standard particle swarm optimization (PSO) and a social
learning PSO for selected benchmark functions and demonstrated its
competitiveness when compared with the state-of-the-art algorithms under a
limited computational budget. Li et. al [50] presented a surrogate-assisted
particle swarm optimization algorithm for computationally expensive problems
where two criterion’s were applied in tandem to select candidates for exact
evaluations. They include a performance-based criterion used to exploit the
current global best and accelerate the convergence rate and a distance-based
uncertainty-based criterion is used to enhance exploration that does not
consider the fitness landscape of different problems. The results demonstrated
better performance over several state-of-the-art algorithms for benchmark
functions and propeller design problem was also used. Chen et. al [51]
presented hierarchical surrogate-assisted differential evolution algorithm for
high-dimensional expensive optimization problems with radial basis function
network and benchmark functions were used and further application to oil
reservoir production optimization problem provided promising results. Yi et.
al [52] presented an on-line variable-fidelity surrogate-assisted harmony
search algorithm with multi-level screening strategy that outperformed related
methods and showed promising performance for expensive engineering design
optimization that featured design for a long cylindrical gas-pressure vessel.

Further work using swarm optimisation is part of latest research. Li et. al
[53] presented ensemble of surrogates assisted particle swarm optimization of
medium scale expensive problems which used multiple trial positions for each
particle and selected the promising positions by using the superiority and
uncertainty of the ensemble simultaneously. In order to feature faster
convergence and to avoid wrong global attraction of models, the optima of two
surrogates that featured polynomial regression and radial basis function
models were evaluated in the convergence state of particles. Liao et. al [54]
presented multi-surrogate multi-tasking optimization of expensive problems to
accelerate the convergence by regarding the two surrogates as two related
tasks. Therefore, two optimal solutions found by the multi-tasking algorithm
was evaluated using the real expensive objective function, and both the global
and local models are updated until the allowed computational budget is
exhausted. The results indicated competitive performance with faster
convergence that scales well with the increase in problem dimension for
solving computationally expensive single-objective optimization problems. Dong
et. al [55] presented surrogate-assisted grey wolf optimization for high-
dimensional and computationally expensive black-box problems that featured
RBF-assisted meta-heuristic exploration. The RBF featured knowledge mining
that includes a global search carried out using the grey wolf optimization and
a local search strategy combining global and multi-start local exploration.
The method obtained superior computation efficiency and robustness
demonstrated by comparison tests with benchmark functions. Chen et. al [51]
presented efficient hierarchical surrogate-assisted differential evolution for
high-dimensional expensive optimization using global and local surrogate model
featuring RBF network with an application to an oil reservoir production
optimization problem. The results show that the method was effective for most
benchmark functions and gave promising performance for reservoir production
optimization problem.

Some of the prominent examples for surrogate assisted approaches in the Earth
sciences include modeling water resources [48, 56], computational oceanography
[57], atmospheric general circulation models [58], carbon-dioxide storage and
oil recovery [59], and debris flow models [60].

###  2.3 Landscape evolution models and Bayeslands

Landscape evolution models (LEMs) use different climate and geophysical
aspects such as tectonics or climate variability [61, 62, 63, 64, 65] and
combine empirical data and conceptual methods into a set of mathematical
equations that form the basis for driving model simulation. Badlands (basin
and landscape dynamics) [63, 66] is a LEM that simulates landscape evolution
and sediment transport/deposition [67, 68] with an initial topography exposed
to climate and geological factors over time with given conditions (parameters)
such as the precipitation rate and rock erodibility coefficient. The major
challenge is in estimating the climate and geological parameters and those
that are linked with the initial topography since they can range millions of
years in geological timescale depending on the problem. A way is to develop an
optimisation framework that utilize limited data. Since gradient information
is not available, Badlands can be seen as a black-box optimisation model where
the unknown parameters need to be found via optimisation or Bayesian
inference. So far, the problem has been approached with Bayesian inference
that employs MCMC sampling for estimation of these parameters in a framework
known as Bayeslands [69].

The Bayeslands framework had limitations due to computational complexity of
Badldands model, hence it was extended using parallel tempering MCMC [70] that
featured parallel computing to enhance computational efficiency and sampling.
Although we used parallel computing with small-scale synthetic Badlands model,
the procedure remained computationally challenging since thousands of samples
were drawn and evaluated. Running a single large scale real-world Badlands
model can take several minutes to hours, and even several days depending on
the area covered by the landscape evolution considered and the span of
geological time considered in terms of millions of years. Therefore, parallel
tempering Bayeslands was further enhanced through surrogate-assisted
estimation. We used developed surrogate-assisted parallel tempering has for
landscape evolution models where a global-local surrogate framework utilised
surrogate training in the main process that managed MCMC replicas running in
parallel [32]. We obtained promising results where prediction performance was
maintained while lowering computational time using surrogates for synthetic
examples that also featured landscape evolution of Tasmania over a million
years.

##  3 Surrogate-assisted distributed swarm optimisation

###  3.1 Particle Swarm Optimisation

Particle swarm optimisation (PSO) is a population based metaheuristic that
improves the population over iterations with a given measure of accuracy known
as fitness [71]. The population of candidate solution in PSO is known as swarm
while the candidate solutions are known as particles that get updated
according to the particle’s position and velocity. In the swarm, each
particle’s movement is typically influenced by its local best known position,
but also guided toward the best known positions that get updated when better
positions are discovered by other particles. Equations 1 and 2 show the
velocity and position update of the particle in a swarm, respectively.

$\mathbf{v}_{t+1}=\alpha\mathbf{v}_{t}+c_{1}\gamma_{1}(\mathbf{x}_{pbest}-\mathbf{x}_{t})+c_{2}\gamma_{2}(\mathbf{x}_{gbest}-\mathbf{x}_{t})$
(1) $\mathbf{x}_{t+1}=\mathbf{x}_{t}+\mathbf{v}_{t+1}$ (2)

where $\boldsymbol{v_{t}}$ $\boldsymbol{x_{t}}$ represent the velocity and
position of a particle at time step $t$, respectively. $c_{1}$ and $c_{2}$ are
the user-defined cognitive and social acceleration coefficients, respectively.
$\gamma_{1}$ and $\gamma_{2}$ are random numbers drawn from uniform $U[0,1]$
distribution, and $\alpha$ is the user-defined inertia weight. There are
several variants in the way the particles get updated which have their own
strengths and limitations for different types of problems [71, 72, 73, 74,
75].

###  3.2 Surrogate assisted Distributed framework

We use a canonical particle update method [71] where the swarms are executed
in parallel computing framework as a separate process.

####  3.2.1 Inter Process Communication

In our distributed swarm framework that used parallel computing, we exchange
selected swarm particles after certain number of generations with inter-
process communication. During the exchange, we replace 20 percent of the
weaker particles by stronger ones from other swarm. This provides another mode
of search space exploration and exploitation to our optimization problem.
Afterwards, the process continues where local swarms create particles with new
position and velocity as shown in Figure 1 . Hence, we feature distributed
swarms and parallel processing for better diversity and computational
complexity. We execute distinct parallel process for respective swarms with
central processing units (CPUs) if enough are present or be shared. Each
processing core features the optimisation function which is synthetic
benchmark problem made to be computationally expensive or application problem
with data and swarms which communicates with the master or manager process.

####  3.2.2 Use of Surrogate Model

Suppose that the true function or model is represented as $F=g(x)$ where $g()$
is the function and $x$ is a solution or particle from a given swarm. Our
surrogate model outputs pseudo-fitness $\hat{F}=\hat{g}(x)$ that would give an
approximation of the true function via $F=\hat{F}+e$; where $e$ represents the
difference between the surrogate and true function. The surrogate model gives
an estimate using the pseudo-fitness for replacing true-function when required
by the framework. $S_{prob}$ is an important hyperparameter which controls use
of surrogate in the prediction, neither want it to be too high in case we are
not very confident about our surrogate model as then the optimization will
become random nor we want it to be too low as then the process will become
time consuming. Hence we need to tune this parameter. Surrogate model is
trained accumulating the data from all the swarms, i.e. the evaluated input
$\bf{x}_{i,s}$ and associated true-fitness $F_{i,s}$ pairs; where $s$
represents the particle and $i$ represents the swarm. To take benifit from the
surrogate in the optimization process its very crucial to manage the surrogate
training and surrogate use. We cannot use surrogate from the very start of the
optimization as its predictions will be random, neither we can wait too long
as the purpose of computational efficiency will be defeated then. Hence to
manage the surrogate training and its use, an important hyperparameter $\psi$
is used, $\psi$ is the surrogate interval measured in terms of number of
generations, it is the interval after which all the collected data is used to
train the model, this updated model is used as the surrogate till the next
interval is reached. The collected input features ($\Phi$) combined with the
true fitness $\lambda$ create $\theta$ for the surrogate model.

$\displaystyle\Phi$ $\displaystyle=$
$\displaystyle([\bf{x}_{1,s},\ldots,\bf{x}_{1,s+\psi}],\ldots,[\bf{x}_{M,s},\ldots,\bf{x}_{M,s+\psi}])$
$\displaystyle\lambda$ $\displaystyle=$
$\displaystyle([F_{1,s},\ldots,F_{1,s+\psi}],\ldots,[F_{M,s},\ldots,F_{M,s+\psi}])$
$\displaystyle\Theta$ $\displaystyle=$ $\displaystyle[\Phi,\lambda]$ (3)

where $\bf{x}_{i,s}$ represents the given particle from the swarm, $s$,
$F_{i,s}$ is the output from the true-fitness, and $M$ is the number of
swarms. Therefore, the surrogate training dataset ($\Theta=[\Phi,\lambda]$) is
made up of input features ($\Phi$) and response ($\lambda$) for the particles
that get collected in each surrogate interval ($s+\psi$). Our pseudo-fitness
is given by $\hat{y}=\hat{f}(\Theta)$.

####  3.2.3 Proposed Framework

In Algorithm 1, we present further details about the methodology that features
surrogate assisted optimisation using distributed swarms.

We implement the algorithm using distributed computation over CPU cores, as
shown in Algorithm 1, the manager process is in blue where inter-process
communication among swarms takes place which exchange parameters at regular
intervals (given by $\phi$) and surrogate model is trained at regular
intervals ($\psi$). The parallel swarms have been highlighted in pink in
Algorithm 1 Stage 0 features the initialization of particles in the swarm. We
begin optimisation process by initialising all the respective swarms (Stage
0.1) in the ensemble with random real numbers in a range depending on the
optimisation function or model. Once the swarms are initialised, we begin the
evolution (optimisation) by first evaluating the particles in the swarms using
the fitness (objective) function. We note that we update the best particle and
best fitness for each of the respective swarms in the ensemble afterwards.
Once these basic operations are done, we begin the evolution process where we
create a new set of swarms for the next generation by velocity and position
update (Stage 1.1).

The crux of the method is when we consider whether to evaluate the fitness
function (true fitness) or to use the surrogate model of the fitness function
(pseudo-fitness) when computing fitness values. Stage 1.2 shows how to update
the fitness using either surrogate or true fitness of particle depending on
the interval and $S_{prob}$, initially till the very first surrogate interval
is not reached all evaluations will be from true function. In Stage 1.3, we
calculate the moving average of past three fitness values for a particular
particle by $F_{past}$ = mean($F_{g-1},F_{g-2},F_{g-3}$) to combine with
surrogate model prediction (Stage 1.4). In Stage 1.5 and 1.6, we calculate
actual fitness and save the values for future surrogate training . Our swarm
particle update depends on $x_{pbest}^{i}$ and $x_{gbest}^{i}$ , it may happen
due to poor surrogate performance some bad particles get high fitness score In
order to avoid this issue, we ensure that $x_{pbest}^{i}$ and $x_{gbest}^{i}$
are from true fitness evaluation. In Stage 2.0, given a regular interval
($\phi$), we prepare exchange of selected particles with neighboring swarms
where we replace a given percentage of weak particles (given by fitness
values) to ensure elitism.

We note that before we consider the use of pseudo-fitness, we need to train
the surrogate model with same training data which is created from the true
fitness. Hence we need to collect the training data for the surrogate model
from all the swarms in the ensemble. Hence, in Stage 3.0, the algorithm uses
surrogate training data collected from Stage 1.6 ( $\Theta$ as shown in
Equation 3). In Stage 4, the algorithm trains the surrogate model in the
manager process with data from Stage 3. The knowledge from trained surrogate
model is then used in the fitness estimation as shown in Stage 1.4.

Stage 5.0 executed the termination condition where the algorithm signals the
manager process to decrement number of swarms alive to terminate the swarm
process, when maximum number of evaluations has been reached ($T_{max}$). We
use neural network based surrogate model with Adam optimisation [65]. In order
to validate the performance of the algorithm, we measure the quality of the
surrogate estimate using the root mean squared error (RMSE):

$RMSE=\sqrt{\frac{1}{N}\sum_{i=1}^{N}(F_{i}-\hat{F_{i}})^{2}}$

where $F_{i}$ and $\hat{F_{i}}$ are the true and pseudo-fitness values,
respectively. $N$ denotes the number of times the surrogate model is used for
estimation.

Figure 1 provides visual description of the proposed algorithm. Multiple
swarms consisting of updating population are managed using the manager process
which keeps control of particle exchange and surrogate model update. All the
details are marked in the figure.

Result: Fitness Score

* Set the number of swarms ($M$) in ensemble as $alive$; $alive=M$

* Define the swarm population size ($pop\\_size$), limits of parameters ($maxx,minx$), number of swarm processes ($M$), surrogate interval ($\psi$), swap-interval ($\phi$), and maximum swarm evaluations ($T_{max}$).

while _( $alive\neq 0$_ do

Prepare manager process to execute swarms in parallel cores

for _each $m$ until $M$_ do

$evals=0$

Stage 0: Initialization of parameters

for _each $i$ until $pop\\_size$_ do

0.1: Initialize the position and velocity of particle $x_{i}$ , $v_{i}$
respectively given the bound.

0.2: Initialize $x_{pbest}^{i}=x_{i}$, $\delta_{i}=F(x_{i})$,
$\delta_{best}^{i}=\delta_{i}$ .

end Compute $\delta_{s}$, $x_{gbest}$ among the swarm particles with best
fitness and its corresponding parameter.

while _( $evals<T_{max}$)_ do

Stage 1.0: PSO

for _each $v$ until $\psi$_ do

for _each $k$ until $\phi$_ do

for _each $i$ until $pop\\_size$_ do

1.1 Update $x_{i}$ , $v_{i}$ via PSO update rule as in equation 1 and 2

Check the parameters bounds using $minx$ and $maxx$.

1.2 Estimate $F_{i}$:

Draw $\kappa$ from a Uniform distribution [0,1]

if _$\kappa <$ $S_{prob}$ and $evals>$ $\psi$ _ then

Predict $F_{pseudo}$ using surrogate model to get $F_{i}$

1.3 $F_{past}$ = mean of fitness values of particle $i$ for last three
generations

1.4 $F_{i}$ = 0.5 * $F_{pseudo}$ \+ 0.5 * $F_{past}$

else

1.5 $F_{true}$ = true-fitness for the particle (F($x_{i}$))

1.6 Save $F_{i}$ = $F_{true}$

end

1.6.1 For particle $i$, $F_{local}=F_{i}$

if _$F_{i} <\delta_{best}^{i}$_ then

if _$F_{i}$ is evaluated using surrogate model_ then

$actual_{err}=F(x_{i})$ ; ($TEE$ += 1)

Now $if(actual_{err}<\delta_{best}^{i})$ : $\delta_{best}^{i}\leftarrow
actual_{err}$, $x_{pbest}^{i}\leftarrow x_{i}$

else

$\delta_{best}^{i}\leftarrow F_{i}$ , $x_{pbest}^{i}\leftarrow x_{i}$

end

if _$\delta_{best}^{i} <\delta_{s}$_ then

$\delta_{s}\leftarrow\delta_{best}^{i}$ , $x_{gbest}\leftarrow x_{i}$

end

end

Increment $evals+=pop\\_size$

end

end

Stage 2.0: Neighbouring swarm exchange:

2.1 Use a swapping probability $\beta$ and draw $b$ from a Uniform
distribution [0,1]

if _$b\leq\beta$_ then

2.2 Signal() manager process

2.3 Exchange some parameters with neighbouring swarms,
$\theta_{m}\leftrightarrow\theta_{m+1}$

end

end Stage 3.0: Collect data $\Theta_{m}$ which features history of parameters
$\Phi$($x$) and fitness response $\lambda$ ( $F_{i}$ ) using data collected
from Stage 1.6

Stage 4.0: Train global surrogate model

for _each island_ do

4.1 Get swarm training data $\Theta_{m}$ from Stage 3.0

end 4.2 Train surrogate model,
$\Upsilon=[\Theta_{1},\Theta_{2},...\Theta_{M}]$ .

4.3 Store surrogate model, $\Psi$

end

Stage 5.0: Signal()the manager process

5.1 Decrement the number of swarms $alive$

end

end Stage 6: Combine predictions for particle in the swarm with best fitness
from different islands in the ensemble.

Algorithm 1 Surrogate assisted swarm optimization for computationally
expensive models.

Figure 1: Surrogate-assisted distributed swarm optimisation features
surrogates to estimate the fitness of expensive models or functions.

##  4 Application: Landscape evolution models

---

(a) Continental-Margin

Figure 2: Location of (a) Continental-Margin problem shown taken from South
Island of New Zealand

Other than optimizing mathematical functions many a time we are required to
evaluate a score using some simulation or other computationally expensive
process, an example of which are Landscape evolution models (LEMs). They
typically require an initial topography that is evolved over given time with
given geological and climate conditions such as rock erodibility and
precipitation [63]. LEMs are used to model and understand the landscape and
basin evolution back in time over millions of years showing surface process
such as formation of river systems and erosion/deposition where there is
movement of sediments from source (mountains) to sink (basins) [29]. LEMs help
geologists and paleoclimate scientists understand the evolution of the planet
and climate history over millions to billions of years; however, there are
major challenges when it comes to data. LEMs generally require data regarding
paleoclimate processes which is unavailable and hence we need to estimate them
with methods such as Bayesian inference [76, 77]. There has not been any work
in the literature where optimisation methods can be used to estimate the
unknown parameters for LEMs, which is the focus of this study. We note that
typically, LEMs are computationally very expensive which is dependent on the
resolution of the study area (points/kilometer) and how how the evolution goes
back in time (millions of years). Hence, large scale study areas can take from
hours to days to run a single model even while using parallel computing [63].
There is typically no gradient information in case of Badlands LEM used for
this study, and hence estimating the parameters is a challenge.

In order to demonstrate the optimisation procedure, we select problems where
synthetic initial topography is created using present day topography or
synthetically generated and used in our previous work [76, 77]. The selected
LEM features a continental margin (CM) problem that is selected taking into
account computational time of a single model run as it takes less than three
seconds on a single central processing unit (CPU). The CM problem initial
topography is selected from present-day South Island of New Zealand as shown
in Figure 2 which covers 136 by 123 kilometers. We provide a visualization of
the initial and final topographies along with erosion/deposition map for CM
problem. in Figure 3. The CM features six free parameters (Table 1). Notable
feature in all three problems is that they model both the elevation and
erosion/deposition topography. We use the initial topography (Figure 3) and
true values in in Table 1 and run the Badlands LEM simulating 1 million years
to synthetically generate ground-truth topography. We then create a fitness
function with the ground-truth topography and set experiments so that the
proposed optimisation methods can get back the true values. The details about
the fitness function is given in the following section.

---

(a) CM initial topography

(b) CM ground-truth topography

(c) CM erosion/deposition map

Figure 3:  The figures shows the initial and evolved ground-truth topography and erosion/deposition after one million years for CM problems taken from [76, 77]. The x-axis denotes the latitude and y-axis denotes the longitude from location given in 2, and the elevation is given in meters which is further shown as a colour-bar.  Issue | Rainfall (m/a) | Erod. | n-value | m-value | c-marine | c-surface | Uplift (mm/a)
---|---|---|---|---|---|---|---
True-values | 1.5 | 5.0-e06 | 1.0 | 0.5 | 0.5 | 0.8 | -
Limits | [0,3.0 ] | [3.0-e06, 7.0-e06] | [0, 2.0] | [0, 2.0] | [0.3, 0.7] | [0.6, 1.0] | -
Table 1: True values and limits of parameters.

###  4.1 Fitness function

The Badlands LEM produces a simulation of successive time-dependent
topographies; however, only the final topography $\mathbf{D}_{T}$ is used for
topography fitness since no successive ground-truth data is available. The
sedimentation (erosion/deposition) data is typically used to ground-truth the
time-dependent evolution of surface process models that include sediment
transportation and deposition [63, 69].

We adapt the fitness function from the likelihood function used in our
previous work that used Bayesian inference via Markov Chain Monte-Carlo
sampling algorithm for parameter estimation in Badlands model [76]. The
initial topography is given as $\boldsymbol{D}_{0}$ with
$\boldsymbol{D}_{0}=(D_{0,s_{1}}\ldots,D_{0,s_{n}})$, where $s_{i}$
corresponds to the location which is given by the latitude $u_{i}$ and
longitude $v_{i}$ (eg. Figure 2). Hence, our topography fitness function
$F_{topo}$ for the topography at final time $t=T$ is given by

$F_{topo}(\boldsymbol{\theta})=\sqrt{1/N\sum_{i=1}^{N}\left(D_{s_{i},T}-f_{s_{i},T}(\boldsymbol{\theta})\right)^{2}}$
(4)

where $\nu$ is the number of observations.

The sediment erosion/deposition values at time ($\mathbf{z}_{t}$) are
simulated (predicted) by the Badlands model when given set of parameters,
$\boldsymbol{\theta}$. The sediment fitness $F_{sed}(\boldsymbol{\theta})$ is
given below

$F_{sed}(\boldsymbol{\theta})=\sqrt{1/(T+J)\sum_{t=1}^{T}\sum_{j=1}^{J}\left(z_{s_{j},t}-g_{s_{j},t}(\boldsymbol{\theta})\right)^{2}}$
(5)

##  5 Experiments and Results

In this section, we present analysis of the proposed methodology first on
synthetic benchmark functions then on Badlands(LEM). The experiments consider
a wide range of performance measures which includes optimization performance
in terms of fitness score, computational time saved by the surrogate model,
accuracy of the surrogate model.

###  5.1 Experiment design

We provide the experimental design and parameter setting for our experiment as
follows. We implement distributed swarm optimization using parallel computing
and inter-process communication where the swarms can have separate process and
exchange solutions (particles) with Python multi-processing library
111https://github.com/sydney-machine-learning/surrogate-assisted-distributed-
swarms.

The synthetic benchmark optimisation functions are given in Equation 6
(Spherical), Equation 7 (Ackley), Equation 8 (Rastrigin), and Equation 9
(Rosenbrock). Note that these functions are chosen due to different level of
difficulty in optimisation and difference in nature of their error or fitness
landscape. Spherical function is considered to be an relatively easier
optimisation problem since it does not have interacting variables. Ackley and
Rastrigin are known to be many local minima functions, while Rosenbrock is
known as the valley shaped function.

${f(\textbf{x})=f(x_{1},x_{2},...,x_{n})={\sum_{i=1}^{n}x_{i}^{2}}}$ (6)

$f(\textbf{x})=\sum_{i=1}^{n}[b(x_{i+1}-x_{i}^{2})^{2}+(a-x_{i})^{2}$ (7)

${f(\textbf{x})=-a*exp(-b\sqrt{\frac{1}{n}\sum_{i=1}^{n}x_{i}^{2}})-exp(\frac{1}{n}\sum_{i=1}^{n}cos(cx_{i}))+a+exp(1)}$
(8)

${f(\textbf{x})=10n+\sum_{i=1}^{n}(x_{i}^{2}-10cos(2\pi x_{i}))}$ (9)

We use (M=8) swarms which run as parallel processes (swarms) that inter-
communicate with each other at regular interval where we exchange of a subset
of particles (best 20% from the population) with the neighbouring island
(worst 20% of the population). We use the population size of 20 particles per
swarm,the inertia weight of 0.729 with social and cognitive
coefficients(c1=1.4 and c2=1.4). We determined these values for parameters in
the trial experiments by taking into account the performance on different
problems with a fixed number of evaluations. We used minimum and maximum bound
on the parameters as described in the table, respectively. We ran experiments
for 30 dimensions (30D) and 50D versions of the respective benchmark problems
where the total number of function evaluations were set to 100 000 and 200
000, respectively. In the case of Badlands model, we used 10 000 function
evaluations for our problem instance, CM.

We used the PyTorch machine learning library 222PyTorch: https://pytorch.org/
for implementing the surrogate model that uses a neural network model with
Adam learning. The surrogate neural network model architecture is given in
Table 2 where the input dimensions are defined by the dimension of the problem
used i.e. 30D or 50D for two different instances of the synthetic benchmark
problems. In a similar way, we extend our approach for optimization in
Badlands model where we have CM problems, featuring 6D search space. The first
and second hidden layers of the of the neural network based surrogate model is
shown in Table Table 2. We note that the output layer in the model for all
problems contains only one unit which provides the estimated fitness score.

Table 2:  Model Architecture for surrogate model Problem |  Hidden($h_{1}$) |  Hidden($h_{2}$) |
---|---|---|---
Rosenbrock 30D | 30 | 15 |
Rastrigin 30D | 30 | 15 |
Ackley 30D | 30 | 15 |
Spherical 30D | 30 | 15 |
Rosenbrock 50D | 50 | 25 |
Rastrigin 50D | 50 | 25 |
Ackley 50D | 50 | 25 |
Spherical 50D | 50 | 25 |
Badlands | 20 | 10 |

###  5.2 Results for synthetic benchmark functions

We first present the results (fitness score) for different problems using
canonical serial PSO, distributed PSO (D-PSO) and surrogate-assisted
distributed swarm optimisation (SD-PSO) as shown in Table 3. We use a fixed
surrogate probability ($s_{prob}=0.5$)and present results featuring mean,
standard deviation (std), best and worst performance for over 30 independent
experimental runs with different random initialisation in swarms as shown in
Table 3. Note the lower fitness scores provides better performance.

We find significant reduction in elapsed time for all the problems in Table 3.
We note that a 0.05 seconds time delay was added to the respective problems to
make them slightly computationally expensive to depict real-world problems
that could take more time. When we consider serial PSO with D-PSO and SD-PSO
in the terms of optimisation performance given by the fitness score, we find
that D-PSO improves the PSO significantly for 30D and 50D cases of Rosenbrock,
Spherical and Rastrigin problems. We also see improvement in Ackley problem,
but its not as large as in previous problems. We find the variation in the
results given by the standard deviation is lowered highly with D-PSO which
shows it is more robust to initialization and has the ability to provide more
definitive solution.

Moreover, we compare the results of the SD-PSO with D-PSO and observe that
overall we get similar or better results with SD-PSO, with reduction of
computation time due to the use of surrogates. Despite the use of surrogate-
based fitness estimation, we observe that the fitness-score has not greatly
depreciated. It is interesting to note in Ackley 30 and 50D case, addition of
surrogates is improved the fitness score. We see a 30% reduction in
computational time , it is expected that this reduction will increase if we
had more time delay (rather than 0.05 seconds) which will will see in
experiments to follow.

The RMSE of prediction of fitness by surrogate model is shown in Table 4. We
notice large RMSE values for Rosenbrock problem when compared to the rest. In
surrogate prediction accuracy (RMSE) given in Figure 5, we observe a constant
reduction in RMSE with intervals (along x-axis) in Ackley (30D and 50D), and
Spherical (30D) model functions. In other problems, the RMSE is lower towards
the end, but trend is not that smooth. We note that in Rosenbrock model, there
exists an interval where our surrogate performs poorly causing major decrease
in accuracy.

Figure 6 and 7 provide a visual description surrogates prediction quality in
the evolution process. We show the bar plots for mean values with 95%
confidence interval (shown by error bars) of actual fitness and pseudo-fitness
at regular intervals for different benchmark problems. We notice that
Rastrigin and Ackley problems have a better surrogate prediction with better
confidence intervals when compared to the Spherical problem. Finally, we
evaluate the effect of surrogate probability for Rosenbrock and Rastrigin 30D
problems. Figure 4 (Panel a) provides graphical analysis of how our fitness
score change with varying surrogate probabilities and computational time is
shown as well (Panel b). We observed that the computational time decreases
linearly with surrogate probability. On the other hand, the fitness score
degrades with surrogate probability, but with an elbow shaped curve and a
trade off can exist between time and optimization performance, which is at
surrogate probability of 0.5.

---
Figure 4:  The figure shows effect of surrogate probability on fitness score and computational time  Table 3: Optimization Results on benchmark functions Problem | Method | Fitness Score | Elapsed Time
---|---|---|---
|  | [mean std best worst] | (minutes)
Rosenbrock 30D | PSO | 183.95 109.94 60.88 455.01 | 84.39
| D-PSO | 108.78 38.0 55.0 178.0 | 11.22
| SD-PSO (0.5) | 149.7 50.91 67.92 270.18 | 7.113
Rosenbrock 50D | PSO | 1285.41 559.54 432.89 2519.02 | 84.37
| D-PSO | 1269.03 390.0 78.0 2084.0 | 11.24
| SD-PSO(0.50) | 2449.32 1093.8 685.48 5298.79 | 7.114
Spherical 30D | PSO | 9.06 9.22 1.43 37.85 | 84.35
| D-PSO | 4.79 3.0 1.0 12.0 | 11.31
| SD-PSO | 1.82 0.74 0.9 3.33 | 7.112
Spherical 50D | PSO | 252.96 102.69 67.75 426.74 | 84.38
| D-PSO | 175.2 45.0 118.0 256.0 | 11.23
| SD-PSO(0.50) | 111.85 48.36 34.05 239.828 | 7.113
Rastrigin 30D | PSO | 72.33 19.72 27.05 112.21 | 84.38
| D-PSO | 62.61 12.0 41.0 90.0 | 11.20
| SD-PSO(0.50) | 77.23 12.57 52.96 107.26 | 7.113
Rastrigin 50D | PSO | 188.34 33.83 130.96 257.53 | 84.39
| D-PSO | 180.7 24.0 122.0 215.0 | 11.1
| SD-PSO(0.50) | 243.73 25.1 192.2 307.82 | 7.114
Ackley 30D | PSO | 2.32 0.39 1.66 3.03 | 84.37
| D-PSO | 1.83 0.0 1.0 2.0 | 11.21
| SD-PSO(0.50) | 1.66 0.33 0.99 2.24 | 7.112
Ackley 50D | PSO | 3.8 0.38 3.24 4.58 | 84.38
| D-PSO | 3.42 0.0 3.0 4.0 | 11.21
| SD-PSO(0.50) | 3.34 0.23 2.76 3.73 | 7.113

---

(a) Ackley

(b) Spherical

(c) Rastrigin

(d) Rosenbrock

(e) Badlands model

Figure 5:  RMSE of surrogate training during different surrogate intervals.  Table 4: Surrogate Accuracy  Problem | Surrogate Prediction RMSE | Surrogate Training RMSE
---|---|---
|  | [mean std]
Rosenbrock 30D | 368.33 | 7.68e+02 1.97e+03
Rosenbrock 50D | 284.73 | 1.09e+03 3.45e+03
Spherical 30D | 0.44 | 1.20e-01 2e-01
Spherical 50D | 0.61 | 7.10e-01 8.75e-01
Rastrigin 30D | 0.46 | 2.88e-01 2.54e-01
Rastrigin 50D | 0.54 | 7.33e-01 1.35e+00
Ackley 30D | 0.08 | 1.15e-02 4.53e-03
Ackley 50D | 0.08 | 1.25e-02 4.17e-03
Badlands (CM) | 0.0218 | 2.59e-3 1.95e-3

---

(a) Rosenbrock 30D

(b) Rosenbrock 50D

(c) Spherical 30D

(d) Spherical 50

Figure 6:  Surrogate fitness versus actual fitness score over number of
fitness evaluations (surrogate evaluations) for different problems.

---

(a) Rastrigin 30D

(b) Rastrigin 50D

(c) Ackley 30D

(d) Ackley 50D

Figure 7:  Surrogate fitness versus actual fitness score over number of
fitness evaluations (surrogate evaluations) for different problems.

###  5.3 Results for the Badlands model

Finally, we present results for the case of the Badlands CM problem which is a
6D problem. The results for CM problem highlighting our methods (PSO, D-PSO
and S-DPSO) when compared to previous approaches (PT-Bayeslands and SAPT-
Bayeslands) are shown in Table 5. The results show the computational time and
prediction performance of the Badlands model in terms of elevation and
erosion/deposition RMSE (using Equation 4 and 5, respectively) given the
optimised parameters. The results show the mean and standard deviation from 30
experimental runs from indpenendt initial positions. We see major reduction in
computational time using DPSO when compared to PSO and find consistent
performance in trems of elevation and erosion RMSE. The experiments used
surrogate model at an interval of 10 generations with probability of 0.5. We
observe that the S-DPSO further improved the performance in terms of
computational time and efficient as shown. The RMSE of the estimation of
fitness function by the surrogate model when compared to the actual Badlands
model is shown in Table 4. We note that the RMSE here cant be compared to the
synthetic fitness functions (eg. Rosenbrock) since the fitness function is
completely different. In synthetic fitness functions, there is no data whereas
in Badlands LEM, we use Badlands prediction and ground-truth topography data
to computer the fitness. Figure 5 (Panel e) shows surrogate training accuracy
(RMSE) for different surrogate intervals where we observe a constant
improvement of performance by surrogate model over time (surrogate intervals).
This implies that the surrogate model is improving as it gathers more data
over time.

In Figure 8, we show the change in CM topology over selected time-slices
simulated by Badlands according to the parameters optimized by S-DPSO. The
elevation RMSE in 5 considers the difference between ground-truth topography
given in Figure 3. We notice that visually the final topography (present day)
in Figure 8 (Panel f) resembles Figure 3 (Panel b). Furthermore, we show in
Figure 9 a cross-section (Panel a) for Badlands predicted elevation vs the
ground-truth elevation for final or present day topography. We also show the
bar-plot (Panel b) of predicted vs ground-truth sediment erosion/deposition at
10 selected locations taken from Figure 3 (Panel c). The cross-section and
bar-plots show that the Badlands prediction well resembled the ground-truth
data,respectively. We obverse that the cross-section (Panel a) uncertainty is
higher for certain locations as highlighted. The high uncertainty is in area
of high slope below seal-level which is reasonable given effect of sediment
flow due to precipitation.

---

/

/

/

Figure 8: This is how topology evolves according to our badlands model optimized using S-DPSO. Note that the distance over x-axis and y-axis is given in kilometers (km) and the elevation is given in meters (m).  Table 5: Badlands performance  Method | Elevation RMSE | Elevation RMSE | Erosion RMSE | Erosion RMSE | Time
---|---|---|---|---|---
| (mean) | (std) | (mean) | (std) | (seconds)
PT-Bayeslands [69], | 70.80 | 10.03 | 35.91 | 11.36 | 3243
SAPT-Bayeslands [70] | 82.0 | 8.23 | 44.33 | 13.37 | 1859
PSO | 60.93 | 3.56 | 18.47 | 4.12 | 18480
D-PSO | 61.61 | 4.46 | 26.97 | 5.39 | 4020
S-DPSO | 61.17 | 4.57 | 23.63 | 5.11 | 3462

---

(a) Cross-section reconstruction by optimizing badlands(CM) using S-DPSO.
Ground Truth is the actual topography whereas Badlands Prediction is the mean
of topography predicted by model optimized using S-DPSO over 30 runs. $95\%$
confidence interval is represented in form of uncertainty.

(b) Sediment prediction using S-DPSO. Ground Truth is the actual topography
whereas Badlands Prediction is the mean of topography predicted by model
optimized using S-DPSO over 30 runs.

Figure 9: Prediction cross-section (Panel a) and sediment erosion/deposition
(Panel b) with uncertainty given as 95 % confidence interval from 30
experimental runs.

##  6 Discussion

Overall, the results show that distributed swarm optimisation can be improved
in performance such as as decrease computation time while retaining
optimisation accuracy (fitness) using our proposed surrogate-assisted
distributed swarm optimisation. We also found that the use of surrogates may
even boost performance in terms of optimisation accuracy as in case of
Spherical, and Ackley functions (Table 3) , and Badlands landscape evolution
model problem (Table 5). This might be due to smooth approximation of the
problem by our surrogate model which prevented the optimization process in
trapping in local minima. We highlight that the Badlands model does not
provide gradient information regarding the parameters and hence only gradient
free optimisation or inference methods can be used. In previous work, we used
MCMC methods (PT-Bayeslands and SAPT-Bayeslands [69, 70]) where the search was
done via random-walk proposal distribution with MCMC replicas running in
parallel. Our results show that the use of meta-heuristic search operators
from particle swarm optimisation provides better search features. The results
motivates the use the proposed methodology for expansive optimisation models,
which can feature other geoscientic models. Further use of surrogate in larger
instances of the Badlands landscape evolution model can provide significant
reduction in computational time.

A major contribution of the methodology has been in the implementation using
parallel computing, which takes into account inter-process communication when
exchanging particles (solutions) during the optimisation process. In our
proposed framework, the surrogate training was implemented in the manager
process and the trained parameters were transferred to the parallel swarm
processes where local surrogate model was used to estimate the fitness of the
particle when required (Figure 1). This implementation had to seamlessly
update the local surrogate model at regular intervals set by the user. We note
that although less than eight parallel swarm processes were used, in large
scale problems, the same implementation can be extended and amended. We note
that in case when the number of parameters in the actual model used for
optimisation significantly increases, different ways of training the surrogate
model can be explored.

Another major contribution from the optimisation process for the case of
landscape evolution model is the inference of the parameters, such as
precipitation values in the Badlands model. Through optimisation, we can
estimate what precipitation values gave rise to the evolution of landscape
which resulted in present day landscape. The landscape evolution model hence
provides a temporal topography map of the geological history of the region
under study, and eventually the planet given data availability which is a
major challenge in geoscience. These topographical maps, along with the
optimised values for geological and climate parameters (such as precipitation
and erodibility) can be very useful to geologies and climate scientists who
study ancient climates.

##  7 Conclusions and Future Work

We presented a surrogate framework that features parallel swarm optimisation
processes and seamlessly integrates surrogate training from the manager
process to enable surrogate fitness estimation. Our results indicate that the
proposed surrogate-assisted optimisation method significantly reduces the
computational time, while retaining solution accuracy. In certain cases, it
also helps in improving the solution accuracy by escaping from local minimum
via the surrogates. Although we used particle swarm optimisation as the
designated algorithm, other optimisation algorithms, such as genetic
algorithms, evolution strategies and differential evolution can also be used.

In future work, the parallel optimisation process could be improved by a
combination of different optimisation algorithms which can provide different
capabilities in terms of exploration and exploitation of the search space. The
proposed framework can also incorporate other benchmark function models,
particularly those that feature constraints and also be applied to discrete
parameter optimisation problems which are expensive computationally.

## Software and Data

We provide an open-source implementation of the proposed algorithm in Python
along with data and sample results 333Surrogate-assisted distributed
evolutionary algo: https://github.com/sydney-machine-learning/surrogate-
assisted-distributed-evo-alg.

## Acknowledgement

The authors would like to thank Ratneel Deo for support during the initial
phase of this research project.

## References

  * [1] L. Davis, Handbook of genetic algorithms (1991).
  * [2] Z. Li, X. Lin, Q. Zhang, H. Liu, Evolution strategies for continuous optimization: A survey of the state-of-the-art, Swarm and Evolutionary Computation 56 (2020) 100694.
  * [3] M. N. Ab Wahab, S. Nefti-Meziani, A. Atyabi, A comprehensive review of swarm optimization algorithms, PloS one 10 (5) (2015) e0122827.
  * [4] W. Fang, J. Sun, Y. Ding, X. Wu, W. Xu, A review of quantum-behaved particle swarm optimization, IETE Technical Review 27 (4) (2010) 336–348.
  * [5] D. Freitas, L. G. Lopes, F. Morgado-Dias, Particle swarm optimisation: A historical review up to the current developments, Entropy 22 (3) (2020) 362.
  * [6] K.-F. Man, K.-S. Tang, S. Kwong, Genetic algorithms: concepts and applications [in engineering design], IEEE transactions on Industrial Electronics 43 (5) (1996) 519–534.
  * [7] A. Ponsich, A. L. Jaimes, C. A. C. Coello, A survey on multiobjective evolutionary algorithms for the solution of the portfolio optimization problem and other finance and economics applications, IEEE Transactions on Evolutionary Computation 17 (3) (2012) 321–344.
  * [8] D. E. Moriarty, A. C. Schultz, J. J. Grefenstette, Evolutionary algorithms for reinforcement learning, Journal of Artificial Intelligence Research 11 (1999) 241–276.
  * [9] P. Larranaga, C. M. H. Kuijpers, R. H. Murga, I. Inza, S. Dizdarevic, Genetic algorithms for the travelling salesman problem: A review of representations and operators, Artificial Intelligence Review 13 (2) (1999) 129–170.
  * [10] P. B. Holden, N. R. Edwards, J. Hensman, R. D. Wilkinson, Abc for climate: dealing with expensive simulators, Handbook of approximate Bayesian computation (2018) 569–95.
  * [11] R. Chandra, D. Azam, A. Kapoor, R. D. Müller, Surrogate-assisted Bayesian inversion for landscape and basin evolution models, Geoscientific Model Development 13 (7) (2020) 2959–2979.
  * [12] Q. Zhang, L. T. Yang, Z. Chen, P. Li, A survey on deep learning for big data, Information Fusion 42 (2018) 146–157.
  * [13] K. M. Bakwad, S. S. Pattnaik, B. Sohi, S. Devi, S. V. Gollapudi, C. V. Sagar, P. Patra, Fast motion estimation using small population-based modified parallel particle swarm optimisation, International Journal of Parallel, Emergent and Distributed Systems 26 (6) (2011) 457–476.
  * [14] S. Sivanandam, P. Visalakshi, Dynamic task scheduling with load balancing using parallel orthogonal particle swarm optimisation, International Journal of Bio-Inspired Computation 1 (4) (2009) 276–286.
  * [15] E. Alba, M. Tomassini, Parallelism and evolutionary algorithms, IEEE transactions on evolutionary computation 6 (5) (2002) 443–462.
  * [16] R. M. Hicks, P. A. Henne, Wing design by numerical optimization, Journal of Aircraft 15 (7) (1978) 407–412.
  * [17] Y. Jin, Surrogate-assisted evolutionary computation: Recent advances and future challenges, Swarm and Evolutionary Computation 1 (2) (2011) 61–70.
  * [18] Y. S. Ong, P. B. Nair, A. J. Keane, Evolutionary optimization of computationally expensive problems via surrogate modeling, AIAA journal 41 (4) (2003) 687–696.
  * [19] Z. Zhou, Y. S. Ong, P. B. Nair, A. J. Keane, K. Y. Lum, Combining global and local surrogate models to accelerate evolutionary optimization, IEEE Transactions on Systems, Man, and Cybernetics, Part C (Applications and Reviews) 37 (1) (2007) 66–76.
  * [20] R. Calandra, A. Seyfarth, J. Peters, M. P. Deisenroth, Bayesian optimization for learning gaits under uncertainty, Annals of Mathematics and Artificial Intelligence 76 (1) (2016) 5–23.
  * [21] Y. S. Ong, P. Nair, A. Keane, K. Wong, Surrogate-assisted evolutionary optimization frameworks for high-fidelity engineering design problems, in: Knowledge Incorporation in Evolutionary Computation, Springer, 2005, pp. 307–331.
  * [22] S. Jeong, M. Murayama, K. Yamamoto, Efficient optimization design method using kriging model, Journal of aircraft 42 (2) (2005) 413–420.
  * [23] A. Samad, K.-Y. Kim, T. Goel, R. T. Haftka, W. Shyy, Multiple surrogate modeling for axial compressor blade shape optimization, Journal of Propulsion and Power 24 (2) (2008) 301–310.
  * [24] S. Greenhill, S. Rana, S. Gupta, P. Vellanki, S. Venkatesh, Bayesian optimization for adaptive experimental design: A review, IEEE Access 8 (2020) 13937–13948.
  * [25] J. Snoek, H. Larochelle, R. P. Adams, Practical bayesian optimization of machine learning algorithms, in: Advances in neural information processing systems, 2012, pp. 2951–2959.
  * [26] B. Shahriari, K. Swersky, Z. Wang, R. P. Adams, N. De Freitas, Taking the human out of the loop: A review of bayesian optimization, Proceedings of the IEEE 104 (1) (2015) 148–175.
  * [27] R. Chandra, R. D. Müller, D. Azam, R. Deo, N. Butterworth, T. Salles, S. Cripps, Multicore parallel tempering Bayeslands for basin and landscape evolution, Geochemistry, Geophysics, Geosystems 20 (11) (2019) 5082–5104.
  * [28] T. Coulthard, Landscape evolution models: a software review, Hydrological processes 15 (1) (2001) 165–173.
  * [29] A. Chen, J. Darbon, J.-M. Morel, Landscape evolution models: A review of their fundamental equations, Geomorphology 219 (2014) 68–86.
  * [30] P. Bishop, Long-term landscape evolution: linking tectonics and surface processes, Earth Surface Processes and Landforms: the Journal of the British Geomorphological Research Group 32 (3) (2007) 329–365.
  * [31] A. Temme, J. Baartman, J. Schoorl, Can uncertain landscape evolution models discriminate between landscape responses to stable and changing future climate? a millennial-scale test, Global and Planetary Change 69 (1-2) (2009) 48–58.
  * [32] R. Chandra, D. Azam, A. Kapoor, R. D. Müller, Surrogate-assisted bayesian inversion for landscape and basin evolution models, Geoscientific Model Development 13 (7) (2020) 2959–2979.
  * [33] A. Lin, W. Sun, H. Yu, G. Wu, H. Tang, Global genetic learning particle swarm optimization with diversity enhancement by ring topology, Swarm and Evolutionary Computation 44 (2019) 571 – 583.
  * [34] H. Liu, X.-W. Zhang, L.-P. Tu, A modified particle swarm optimization using adaptive strategy, Expert Systems with Applications 152 (2020) 113353.
  * [35] M. Roshanzamir, M. A. Balafar, S. N. Razavi, A new hierarchical multi group particle swarm optimization with different task allocations inspired by holonic multi agent systems, Expert Systems with Applications 149 (2020) 113292\.
  * [36] B. S. P. Mishra, S. Dehuri, Parallel computing environments: A review, IETE Technical Review 28 (3) (2011) 240–247.
  * [37] C. Leopold, Parallel and Distributed Computing: A survey of models, paradigms, and approaches, Wiley New York, 2001.
  * [38] Y.-J. Gong, W.-N. Chen, Z.-H. Zhan, J. Zhang, Y. Li, Q. Zhang, J.-J. Li, Distributed evolutionary algorithms and their models: A survey of the state-of-the-art, Applied Soft Computing 34 (2015) 286–300.
  * [39] E. Cantú-Paz, A survey of parallel genetic algorithms, Calculateurs paralleles, reseaux et systems repartis 10 (2) (1998) 141–171.
  * [40] J. M. Hereford, A distributed particle swarm optimization algorithm for swarm robotic applications, in: 2006 IEEE International Conference on Evolutionary Computation, 2006, pp. 1678–1685.
  * [41] H. Abadlia, N. Smairi, K. Ghedira, Particle swarm optimization based on dynamic island model, in: 2017 IEEE 29th International Conference on Tools with Artificial Intelligence (ICTAI), IEEE, 2017, pp. 709–716.
  * [42] A. de Campos, A. T. Pozo, E. P. Duarte, Parallel multi-swarm pso strategies for solving many objective optimization problems, Journal of Parallel and Distributed Computing 126 (2019) 13 – 33.  doi:https://doi.org/10.1016/j.jpdc.2018.11.008.
URL http://www.sciencedirect.com/science/article/pii/S0743731518308554

  * [43] A. Erdeljan, D. Capko, S. Vukmirovic, D. Bojanic, V. Congradac, Distributed pso algorithm for data model partitioning in power distribution systems, Journal of Applied Research and Technology 12 (5) (2014) 947 – 957.
  * [44] D. Lim, Y. Jin, Y.-S. Ong, B. Sendhoff, Generalizing surrogate-assisted evolutionary computation, IEEE Transactions on Evolutionary Computation 14 (3) (2010) 329–355.
  * [45] A. Giunta, L. Watson, A comparison of approximation modeling techniques-polynomial versus interpolating models, in: 7th AIAA/USAF/NASA/ISSMO Symposium on Multidisciplinary Analysis and Optimization, p. 4758.
  * [46] R. Jin, W. Chen, T. W. Simpson, Comparative studies of metamodelling techniques under multiple modelling criteria, Structural and multidisciplinary optimization 23 (1) (2001) 1–13.
  * [47] A. Díaz-Manríquez, G. Toscano, J. H. Barron-Zambrano, E. Tello-Leal, A review of surrogate assisted multiobjective evolutionary algorithms, Computational intelligence and neuroscience 2016 (2016).
  * [48] S. Razavi, B. A. Tolson, D. H. Burn, Review of surrogate modeling in water resources, Water Resources Research 48 (7) (2012).
  * [49] H. Yu, Y. Tan, J. Zeng, C. Sun, Y. Jin, Surrogate-assisted hierarchical particle swarm optimization, Information Sciences 454 (2018) 59–72.
  * [50] F. Li, W. Shen, X. Cai, L. Gao, G. Gary Wang, A fast surrogate-assisted particle swarm optimization algorithm for computationally expensive problems, Applied Soft Computing 92 (2020) 106303.
  * [51] G. Chen, Y. Li, K. Zhang, X. Xue, J. Wang, Q. Luo, C. Yao, J. Yao, Efficient hierarchical surrogate-assisted differential evolution for high-dimensional expensive optimization, Information Sciences 542 (2021) 228 – 246.
  * [52] J. Yi, L. Gao, X. Li, C. A. Shoemaker, C. Lu, An on-line variable-fidelity surrogate-assisted harmony search algorithm with multi-level screening strategy for expensive engineering design optimization, Knowledge-Based Systems 170 (2019) 1 – 19.
  * [53] F. Li, X. Cai, L. Gao, Ensemble of surrogates assisted particle swarm optimization of medium scale expensive problems, Applied Soft Computing 74 (2019) 291 – 305.
  * [54] P. Liao, C. Sun, G. Zhang, Y. Jin, Multi-surrogate multi-tasking optimization of expensive problems, Knowledge-Based Systems 205 (2020) 106262.
  * [55] H. Dong, Z. Dong, Surrogate-assisted grey wolf optimization for high-dimensional, computationally expensive black-box problems, Swarm and Evolutionary Computation 57 (2020) 100713.
  * [56] M. J. Asher, B. F. Croke, A. J. Jakeman, L. J. Peeters, A review of surrogate models and their application to groundwater modeling, Water Resources Research 51 (8) (2015) 5957–5973.
  * [57] R. van der Merwe, T. K. Leen, Z. Lu, S. Frolov, A. M. Baptista, Fast neural network surrogates for very high dimensional physics-based models in computational oceanography, Neural Networks 20 (4) (2007) 462–478.
  * [58] S. Scher, Toward data-driven weather and climate forecasting: Approximating a simple general circulation model with deep learning, Geophysical Research Letters 45 (0) (2018) 1–7.  doi:10.1029/2018GL080704.
URL https://agupubs.onlinelibrary.wiley.com/doi/abs/10.1029/2018GL080704

  * [59] W. Ampomah, R. Balch, R. Will, M. Cather, D. Gunda, Z. Dai, Co-optimization of co2-eor and storage processes under geological uncertainty, Energy Procedia 114 (2017) 6928–6941.
  * [60] M. Navarro, O. P. Le Maître, I. Hoteit, D. L. George, K. T. Mandli, O. M. Knio, Surrogate-based parameter inference in debris flow model, Computational Geosciences (2018) 1–17.
  * [61] K. X. Whipple, G. E. Tucker, Implications of sediment-flux-dependent river incision models for landscape evolution, Journal of Geophysical Research: Solid Earth 107 (B2) (2002) 1–20.
  * [62] G. E. Tucker, G. R. Hancock, Modelling landscape evolution, Earth Surface Processes and Landforms 35 (1) (2010) 28–50.
  * [63] T. Salles, X. Ding, G. Brocard, pybadlands: A framework to simulate sediment transport, landscape dynamics and basin stratigraphic evolution through space and time, PloS one 13 (4) (2018) e0195557.
  * [64] B. Campforts, W. Schwanghart, G. Govers, Accurate simulation of transient landscape evolution by eliminating numerical diffusion: the ttlem 1.0 model, Earth Surface Dynamics 5 (1) (2017) 47–66.
  * [65] J. M. Adams, N. M. Gasparini, D. E. J. Hobley, G. E. Tucker, E. W. H. Hutton, S. S. Nudurupati, E. Istanbulluoglu, The landlab v1.0 overlandflow component: a python tool for computing shallow-water flow across watersheds, Geoscientific Model Development 10 (4) (2017) 1645–1663.
  * [66] T. Salles, L. Hardiman, Badlands: An open-source, flexible and parallel framework to study landscape dynamics, Computers & Geosciences 91 (2016) 77–89.
  * [67] A. D. Howard, W. E. Dietrich, M. A. Seidl, Modeling fluvial erosion on regional to continental scales, Journal of Geophysical Research: Solid Earth 99 (B7) (1994) 13971–13986.
  * [68] D. E. J. Hobley, H. D. Sinclair, S. M. Mudd, P. A. Cowie, Field calibration of sediment flux dependent river incision, Journal of Geophysical Research: Earth Surface 116 (F4) (2011).
  * [69] R. Chandra, D. Azam, R. D. Müller, T. Salles, S. Cripps, Bayeslands: A Bayesian inference approach for parameter uncertainty quantification in badlands, Computers & Geosciences 131 (2019) 89–101.
  * [70] R. Chandra, R. D. Müller, D. Azam, R. Deo, N. Butterworth, T. Salles, S. Cripps, Multicore parallel tempering Bayeslands for basin and landscape evolution, Geochemistry, Geophysics, Geosystems 20 (11) (2019) 5082–5104.
  * [71] J. Kennedy, R. Eberhart, Particle swarm optimization, in: Proceedings of ICNN’95-International Conference on Neural Networks, Vol. 4, IEEE, 1995, pp. 1942–1948.
  * [72] Y. Shi, et al., Particle swarm optimization: developments, applications and resources, in: Proceedings of the 2001 congress on evolutionary computation (IEEE Cat. No. 01TH8546), Vol. 1, IEEE, 2001, pp. 81–86.
  * [73] F. Van den Bergh, A. P. Engelbrecht, A cooperative approach to particle swarm optimization, IEEE transactions on evolutionary computation 8 (3) (2004) 225–239.
  * [74] Z.-H. Zhan, J. Zhang, Y. Li, H. S.-H. Chung, Adaptive particle swarm optimization, IEEE Transactions on Systems, Man, and Cybernetics, Part B (Cybernetics) 39 (6) (2009) 1362–1381.
  * [75] S. Yang, M. Wang, et al., A quantum particle swarm optimization, in: Proceedings of the 2004 Congress on Evolutionary Computation (IEEE Cat. No. 04TH8753), Vol. 1, IEEE, 2004, pp. 320–324.
  * [76] R. Chandra, R. D. Müller, D. Azam, R. Deo, N. Butterworth, T. Salles, S. Cripps, Multi-core parallel tempering bayeslands for basin and landscape evolution, Geochemistry, Geophysics, Geosystems 20 (11) (2019) 5082–5104.  doi:10.1029/2019GC008465.
  * [77] R. Chandra, D. Azam, R. D. Müller, T. Salles, S. Cripps, Bayeslands: A bayesian inference approach for parameter uncertainty quantification in badlands, Computers & Geosciences 131 (2019) 89–101.