Genetic Algorithm

Description

This function allows to use genetic algorithm for numeric global optimization of real-valued functions.

Usage

gena(
  fn,
  gr = NULL,
  lower,
  upper,
  pop.n = 100,
  pop.initial = NULL,
  pop.method = "uniform",
  mating.method = "rank",
  mating.par = NULL,
  mating.self = FALSE,
  crossover.method = "local",
  crossover.par = NULL,
  crossover.prob = 0.8,
  mutation.method = "constant",
  mutation.par = NULL,
  mutation.prob = 0.2,
  mutation.genes.prob = 1/length(lower),
  elite.n = min(10, 2 * round(pop.n/20)),
  elite.duplicates = FALSE,
  hybrid.method = "rank",
  hybrid.par = 2,
  hybrid.prob = 0,
  hybrid.opt.par = NULL,
  hybrid.n = 1,
  constr.method = NULL,
  constr.par = NULL,
  maxiter = 100,
  is.max = TRUE,
  info = TRUE,
  ...
)

Arguments

Details

To find information on particular methods available via pop.method, mating.method, crossover.method, mutation.method, hybrid.method and constr.method arguments please see 'Details' section of gena.population, gena.crossover, gena.mutation, gena.hybrid and gena.constr correspondingly. For example to find information on possible values of mutation.method and mutation.par arguments see description of method and par arguments of gena.mutation function.

It is possible to provide manually implemented functions for population initialization, mating, crossover, mutation and hybridization. For example manual mutation function may be provided through mutation.method argument. It should have the same signature (arguments) as gena.mutation function and return the same object i.e. the matrix of chromosomes of the appropriate size. Manually implemented functions for other operators (crossover, mating and so on) may be provided in a similar way.

By default function does not impose any constraints upon the parameters. If constr.method = "bounds" then lower and upper constraints will be imposed. Lower bounds should be strictly smaller then upper bounds.

Currently the only available termination condition is maxiter. We are going to provide some additional termination conditions during future updates.

Infinite values in lower and upper are substituted with -(.Machine$double.xmax * 0.9) and .Machine$double.xmax * 0.9 correspondingly.

By default if gr is provided then BFGS algorithm will be used inside optim during hybridization. Otherwise Nelder-Mead will be used. Manual values for optim arguments may be provided (as a list) through hybrid.opt.par argument.

Arguments pop.n and elite.n should be even integers and elite.n should be greater then 2. If these arguments are odd integers then they will be coerced to even integers by adding 1. Also pop.n should be greater then elite.n at least by 2.

For more information on the genetic algorithm please see Katoch et. al. (2020).

Value

This function returns an object of class gena that is a list containing the following elements:

• par - chromosome (solution) with the highest fitness (objective function) value.

• value - value of fn at par.

• population - matrix of chromosomes (solutions) of the last iteration of the algorithm.

• counts - a two-element integer vector giving the number of calls to fn and gr respectively.

• is.max - identical to is.max input argument.

• fitness.history - vector which i-th element is fitness of the best chromosome in i-th iteration.

• iter - last iteration number.

References

S. Katoch, S. Chauhan, V. Kumar (2020). A review on genetic algorithm: past, present, and future. Multimedia Tools and Applications, 80, 8091-8126. <doi:10.1007/s11042-020-10139-6>

Examples

## Consider Ackley function

fn <- function(par, a = 20, b = 0.2)
{
  val <- a * exp(-b * sqrt(0.5 * (par[1] ^ 2 + par[2] ^ 2))) +
         exp(0.5 * (cos(2 * pi * par[1]) + cos(2 * pi * par[2]))) -
         exp(1) - a
  return(val)
}

# Maximize this function using classical
# genetic algorithm setup
set.seed(123)
lower <- c(-5, -100)
upper <- c(100, 5)
opt <- gena(fn = fn, 
            lower = lower, upper = upper,
            hybrid.prob = 0,
            a = 20, b = 0.2)
print(opt$par)

# Replicate optimization using hybridization
opt <- gena(fn = fn, 
            lower = lower, upper = upper,
            hybrid.prob = 0.2,
            a = 20, b = 0.2)
print(opt$par)


## Consider Rosenbrock function
fn <- function(par, a = 100)
{
  val <- -(a * (par[2] - par[1] ^ 2) ^ 2 + (1 - par[1]) ^ 2 +
           a * (par[3] - par[2] ^ 2) ^ 2 + (1 - par[2]) ^ 2)
  return(val)
}

# Apply genetic algorithm
lower <- rep(-10, 3)
upper <- rep(10, 3)
set.seed(123)
opt <- gena(fn = fn,
            lower = lower, upper = upper,
            a = 100)
print(opt$par)


# Improve the results by hybridization
opt <- gena(fn = fn,
            lower = lower, upper = upper,
            hybrid.prob = 0.2,
            a = 100)
print(opt$par)


# Provide manually implemented mutation function
# which simply randomly sorts genes. 
# Note that this function should have the same 
# arguments as gena.mutation.
mutation.my <- function(children, lower, upper, 
                        prob, prob.genes, 
                        method, par, iter)
{
  # Get dimensional data
  children.n <- nrow(children)
  genes.n <- ncol(children)
  
  # Select chromosomes that should mutate
  random_values <- runif(children.n, 0, 1)
  mutation_ind <- which(random_values <= prob)
   
  # Mutate chromosomes by randomly sorting
  # their genes
  for (i in mutation_ind)
  {
    children[i, ] <- children[i, sample(1:genes.n)]
  }
  
  # Return mutated chromosomes
  return(children)
}

opt <- gena(fn = fn,
            lower = lower, upper = upper,
            mutation.method = mutation.my,
            a = 100)
print(opt$par)

Constraints

Description

Impose constraints on chromosomes.

Usage

gena.constr(population, method = "bounds", par, iter)

Arguments

Details

If method = "bounds" then chromosomes will be bounded between par$lower and par$upper.

Value

The function returns a list with the following elements:

• population - matrix which rows are chromosomes after constraints have been imposed.

• constr.ind - matrix of logical values which (i, j)-th elements equals TRUE (FALSE otherwise) if j-th jene of i-th chromosome is a subject to constraint.

Examples

# Randomly initialize population
set.seed(123)
population <- gena.population(pop.n = 10,
                              lower = c(-5, -5), 
                              upper = c(5, 5))
                           
# Impose lower and upper bounds constraints
pop.constr <- gena.constr(population, 
                          method = "bounds",
                          par = list(lower = c(-1, 2),
                                     upper = c(1, 5)))
print(pop.constr)

Crossover

Description

Crossover method (algorithm) to be used in the genetic algorithm.

Usage

gena.crossover(
  parents,
  fitness = NULL,
  prob = 0.8,
  method = "local",
  par = NULL,
  iter = NULL
)

Arguments

Details

Denote parents by C^{parent} which i-th row parents[i, ] is a chromosome c_{i}^{parent} i.e. the vector of parameter values of the function being optimized f(.) that is provided via fn argument of gena. The elements of chromosome c_{ij}^{parent} are genes representing parameters values.

Crossover algorithm determines the way parents produce children. During crossover each of randomly selected pairs of parents c_{i}^{parent}, c_{i + 1}^{parent} produce two children c_{i}^{child}, c_{i + 1}^{child}, where i is odd. Each pair of parents is selected with probability prob. If pair of parents have not been selected for crossover then corresponding children and parents are coincide i.e. c_{i}^{child}=c_{i}^{parent} and c_{i+1}^{child}=c_{i+1}^{parent}.

Argument method determines particular crossover algorithm to be applied. Denote by \tau the vector of parameters used by the algorithm. Note that \tau corresponds to par.

If method = "split" then each gene of the first child will be equiprobably picked from the first or from the second parent. So c_{ij}^{child} may be equal to c_{ij}^{parent} or c_{(i+1)j}^{parent} with equal probability. The second child is the reversal of the first one in a sense that if the first child gets particular gene of the first (second) parent then the second child gets this gene from the first (second) parent i.e. if c_{ij}^{child}=c_{ij}^{parent} then c_{(i+1)j}^{child}=c_{(i+1)j}^{parent}; if c_{ij}^{child}=c_{(i+1)j}^{parent} then c_{(i+1)j}^{child}=c_{ij}^{parent}.

If method = "arithmetic" then:

c_{i}^{child}=\tau_{1}c_{i}^{parent}+ \left(1-\tau_{1}\right)c_{i+1}^{parent}

c_{i+1}^{child}=\left(1-\tau_{1}\right)c_{i}^{parent}+ \tau_{1}c_{i+1}^{parent}

where \tau_{1} is par[1]. By default par[1] = 0.5.

If method = "local" then the procedure is the same as for "arithmetic" method but \tau_{1} is a uniform random value between 0 and 1.

If method = "flat" then c_{ij}^{child} is a uniform random number between c_{ij}^{parent} and c_{(i+1)j}^{parent}. Similarly for the second child c_{(i+1)j}^{child}.

For more information on crossover algorithms please see Kora, Yadlapalli (2017).

Value

The function returns a matrix which rows are children.

References

P. Kora, P. Yadlapalli. (2017). Crossover Operators in Genetic Algorithms: A Review. International Journal of Computer Applications, 162 (10), 34-36, <doi:10.5120/ijca2017913370>.

Examples

# Randomly initialize the parents
set.seed(123)
parents.n <- 10
parents <- gena.population(pop.n = parents.n,
                           lower = c(-5, -5), 
                           upper = c(5, 5))
                           
# Perform the crossover
children <- gena.crossover(parents = parents,
                           prob = 0.6,
                           method = "local")
print(children)

Hybridization

Description

Hybridization method (algorithm) to be used in the genetic algorithm.

Usage

gena.hybrid(
  population,
  fitness,
  hybrid.n = 1,
  method,
  par,
  opt.par,
  info = FALSE,
  iter = NULL,
  ...
)

Arguments

Details

This function uses gena.mating function to select hybrids. Therefore method and par arguments will be passed to this function. If some chromosomes selected to become hybrids are duplicated then these duplicates will not be subject to local optimization i.e. the number of hybrids will be decreased by the number of duplicates (actual number of hybrids during some iterations may be lower than hybrid.n).

Currently optim is the only available local optimizer. Therefore opt.par is a list containing parameters that should be passed to optim.

For more information on hybridization please see El-mihoub et. al. (2006).

Value

The function returns a list with the following elements:

• population - matrix which rows are chromosomes including hybrids.

• fitness - vector which i-th element is the fitness of the i-th chromosome.

• hybrids.ind - vector of indexes of chromosomes selected for hybridization.

• counts a two-element integer vector giving the number of calls to fn and gr respectively.

References

T. El-mihoub, A. Hopgood, L. Nolle, B. Alan (2006). Hybrid Genetic Algorithms: A Review. Engineering Letters, 13 (3), 124-137.

Examples

# Consider the following fitness function
fn <- function(x)
{
  val <- x[1] * x[2] - x[1] ^ 2 - x[2] ^ 2
}

# Also let's provide it's gradient (optional)
gr <- function(x)
{
  val <- c(x[2] - 2 * x[1],
           x[1] - 2 * x[2])
}

# Randomly initialize the population
set.seed(123)
n_population <- 10
population <- gena.population(pop.n = n_population,
                              lower = c(-5, -5), 
                              upper = c(5, 5))

# Calculate fitness of each chromosome
fitness <- rep(NA, n_population)
for(i in 1:n_population)
{
  fitness[i] <- fn(population[i, ])
}

# Perform hybridization
hybrids <- gena.hybrid(population = population,
                       fitness = fitness,
                       opt.par = list(fn = fn,
                                      gr = gr,
                                      method = "BFGS",
                                      control = list(fnscale = -1,
                                                     abstol = 1e-10,
                                                     reltol = 1e-10,
                                                     maxit = 1000)),
                       hybrid.n = 2,
                       method = "rank",
                       par = 0.8)
print(hybrids)

Mating

Description

Mating (selection) method (algorithm) to be used in the genetic algorithm.

Usage

gena.mating(
  population,
  fitness,
  parents.n,
  method = "rank",
  par = NULL,
  self = FALSE,
  iter = NULL
)

Arguments

Details

Denote population by C which i-th row population[i, ] is a chromosome c_{i} i.e. the vector of parameter values of the function being optimized f(.) that is provided via fn argument of gena. The elements of chromosome c_{ij} are genes representing parameters values. Argument fitness is a vector of function values at corresponding chromosomes i.e. fitness[i] corresponds to f_{i}=f(c_{i}). Total number of chromosomes in population n_{population} equals to nrow(population).

Mating algorithm determines selection of chromosomes that will become parents. During mating each iteration one of chromosomes become a parent until there are n_{parents} (i.e. parents.n) parents selected. Each chromosome may become a parent multiple times or not become a parent at all.

Denote by c^{s}_{i} the i-th of selected parents. Parents c^{s}_{i} and c^{s}_{i + 1} form a pair that will further produce a child (offspring), where i is odd. If self = FALSE then for each pair of parents (c_{i}^s, c_{i+1}^s) it is insured that c^{s}_{i} \ne c^{s}_{i + 1} except the case when there are several identical chromosomes in population. However self is ignored if method is "tournament", so in this case self-mating is always possible.

Denote by p_{i} the probability of a chromosome to become a parent. Remind that each chromosome may become a parent multiple times. Probability p_{i}\left(f_{i}\right) is a function of fitness f_{i}. Usually this function is non-decreasing so more fitted chromosomes have higher probability of becoming a parent. There is also an intermediate value w_{i} called weight such that:

p_{i}=\frac{w_{i}}{\sum\limits_{j=1}^{n_{population}}w_{j}}

Therefore all weights w_{i} are proportional to corresponding probabilities p_{i} by the same factor (sum of weights).

Argument method determines particular mating algorithm to be applied. Denote by \tau the vector of parameters used by the algorithm. Note that \tau corresponds to par. The algorithm determines a particular form of the w_{i}\left(f_{i}\right) function which in turn determines p_{i}\left(f_{i}\right).

If method = "constant" then all weights and probabilities are equal:

w_{i}=1 => p_{i}=\frac{1}{n_{population}}

If method = "rank" then each chromosome receives a rank r_{i} based on the fitness f_{i} value. So if j-th chromosome is the fittest one and k-th chromosome has the lowest fitness value then r_{j}=n_{population} and r_{k}=1. The relationship between weight w_{i} and rank r_{i} is as follows:

w_{i}=\left(\frac{r_{i}}{n_{population}}\right)^{\tau_{1}}

The greater value of \tau_{1} the greater portion of probability will be delivered to more fitted chromosomes. Default value is \tau_{1} = 0.5 so par = 0.5.

If method = "fitness" then weights are calculated as follows:

w_{i}=\left(f_{i} - \min\left(f_{1},...,f_{n_{population}}\right) + \tau_{1}\right)^{\tau_{2}}

By default \tau_{1}=10 and \tau_{2}=0.5 i.e. par = c(10, 0.5). There is a restriction \tau_{1}\geq0 insuring that expression in brackets is non-negative.

If method = "tournament" then \tau_{1} (i.e. par) chromosomes will be randomly selected with equal probabilities and without replacement. Then the chromosome with the highest fitness (among these selected chromosomes) value will become a parent. It is possible to provide representation of this algorithm via probabilities p_{i} but the formulas are numerically unstable. By default par = min(5, ceiling(parents.n * 0.1)).

Validation and default values assignment for par is performed inside gena function not in gena.mating. It allows to perform validation a single time instead of repeating it each iteration of genetic algorithm.

For more information on mating (selection) algorithms please see Shukla et. al. (2015).

Value

The function returns a list with the following elements:

• parents - matrix which rows are parents. The number of rows of this matrix equals to parents.n while the number of columns is ncol(population).

• fitness - vector which i-th element is the fitness of the i-th parent.

• ind - vector which i-th element is the index of i-th parent in population so $parents[i, ] equals to population[ind[i], ].

References

A. Shukla, H. Pandey, D. Mehrotra (2015). Comparative review of selection techniques in genetic algorithm. 2015 International Conference on Futuristic Trends on Computational Analysis and Knowledge Management (ABLAZE), 515-519, <doi:10.1109/ABLAZE.2015.7154916>.

Examples

# Consider the following fitness function
fn <- function(x)
{
  val <- x[1] * x[2] - x[1] ^ 2 - x[2] ^ 2
}

# Randomly initialize the population
set.seed(123)
pop.nulation <- 10
population <- gena.population(pop.n = pop.nulation,
                              lower = c(-5, -5), 
                              upper = c(5, 5))

# Calculate fitness of each chromosome
fitness <- rep(NA, pop.nulation)
for(i in 1:pop.nulation)
{
  fitness[i] <- fn(population[i, ])
}

# Perform mating to select parents
parents <- gena.mating(population = population,
                       fitness = fitness,
                       parents.n = pop.nulation,
                       method = "rank",
                       par = 0.8)
print(parents)

Mutation

Description

Mutation method (algorithm) to be used in the genetic algorithm.

Usage

gena.mutation(
  children,
  lower,
  upper,
  prob = 0.2,
  prob.genes = 1/nrow(children),
  method = "constant",
  par = 1,
  iter = NULL
)

Arguments

Details

Denote children by C^{child} which i-th row children[i, ] is a chromosome c_{i}^{child} i.e. the vector of parameter values of the function being optimized f(.) that is provided via fn argument of gena. The elements of chromosome c_{ij}^{child} are genes representing parameters values.

Mutation algorithm determines random transformation of children's genes. Each child may be selected for mutation with probability prob. If i-th child is selected for mutation and prob.genes is a vector then j-th gene of this child is transformed with probability prob.genes[j]. If prob.genes is a constant then this probability is the same for all genes.

Argument method determines particular mutation algorithm to be applied. Denote by \tau the vector of parameters used by the algorithm. Note that \tau corresponds to par. Also let's denote by c_{ij}^{mutant} the value of gene c_{ij}^{child} after mutation.

If method = "constant" then c_{ij}^{mutant} is a uniform random variable between lower[j] and upper[j].

If method = "normal" then c_{ij}^{mutant} equals to the sum of c_{ij}^{child} and normal random variable with zero mean and standard deviation par[j]. By default par is identity vector of length ncol(children) so par[j] = 1 for all j.

If method = "percent" then c_{ij}^{mutant} is generated from c_{ij}^{child} by equiprobably increasing or decreasing it by q percent, where q is a uniform random variable between 0 and par[j]. Note that par may also be a constant then all genes have the same maximum possible percentage change. By default par = 20.

For more information on mutation algorithms please see Patil, Bhende (2014).

Value

The function returns a matrix which rows are children (after mutation has been applied to some of them).

References

S. Patil, M. Bhende. (2014). Comparison and Analysis of Different Mutation Strategies to improve the Performance of Genetic Algorithm. International Journal of Computer Science and Information Technologies, 5 (3), 4669-4673.

Examples

# Randomly initialize some children
set.seed(123)
children.n <- 10
children <- gena.population(pop.n = children.n,
                            lower = c(-5, -5), 
                            upper = c(5, 5))
                           
# Perform the mutation
mutants <- gena.mutation(children = children,
                         prob = 0.6,
                         prob.genes = c(0.7, 0.8),
                         par = 30,
                         method = "percent")
print(mutants)

Population

Description

Initialize the population of chromosomes.

Usage

gena.population(pop.n, lower, upper, pop.initial = NULL, method = "uniform")

Arguments

Details

If "method = uniform" then i-th gene of each chromosome is randomly (uniformly) chosen between lower[i] and upper[i] bounds. If "method = normal" then i-th gene is generated from a truncated normal distribution with mean (upper[i] + lower[i]) / 2 and standard deviation (upper[i] - lower[i]) / 6 where lower[i] and upper[i] are lower and upper truncation bounds correspondingly. If "method = hypersphere" then population is simulated uniformly from the hypersphere with center upper - lower and radius sqrt(sum((upper - lower) ^ 2)) via rhypersphere function setting type = "inside".

Value

This function returns a matrix which rows are chromosomes.

References

B. Kazimipour, X. Li, A. Qin (2014). A review of population initialization techniques for evolutionary algorithms. 2014 IEEE Congress on Evolutionary Computation, 2585-2592, <doi:10.1109/CEC.2014.6900618>.

Examples

set.seed(123)
gena.population(pop.n = 10,
                lower = c(-1, -2, -3),
                upper = c(1, 0, -1),
                pop.initial = rbind(c(0, -1, -2),
                                    c(0.1, -1.2, -2.3)),
                method = "normal")

Numeric Differentiation

Description

Numeric estimation of the gradient and Hessian.

Usage

gena.grad(
  fn,
  par,
  eps = sqrt(.Machine$double.eps) * abs(par),
  method = "central-difference",
  fn.args = NULL
)

gena.hessian(
  fn = NULL,
  gr = NULL,
  par,
  eps = sqrt(.Machine$double.eps) * abs(par),
  fn.args = NULL,
  gr.args = NULL
)

Arguments

Details

It is possible to substantially improve numeric Hessian accuracy by using analytical gradient gr. If both fn and gr are provided then only gr will be used. If only fn is provided for gena.hessian then eps will be transformed to sqrt(eps) for numeric stability purposes.

Value

Function gena.grad returns a vector that is a gradient of fn at point par calculated via method numeric differentiation approach using increment eps.

Function gena.hessian returns a matrix that is a Hessian of fn at point par.

Examples

# Consider the following function
fn <- function(par, a = 1, b = 2)
{
  val <- par[1] * par[2] - a * par[1] ^ 2 - b * par[2] ^ 2
}

# Calculate the gradient at point (2, 5) respect to 'par' 
# when 'a = 1' and 'b = 1'
par <- c(2, 5)
fn.args = list(a = 1, b = 1)
gena.grad(fn = fn, par = par, fn.args = fn.args)

# Calculate Hessian at the same point
gena.hessian(fn = fn, par = par, fn.args = fn.args)

# Repeat calculation of the Hessian using analytical gradient
gr <- function(par, a = 1, b = 2)
{
  val <- c(par[2] - 2 * a * par[1],
           par[1] - 2 * b * par[2])
}
gena.hessian(gr = gr, par = par, gr.args = fn.args)

Plot best found fitnesses during genetic algorithm

Description

Plot best found fitnesses during genetic algorithm

Usage

## S3 method for class 'gena'
plot(x, y = NULL, ...)

Arguments

Value

This function does not return anything.

Plot best found fitnesses during genetic algorithm

Description

Plot best found fitnesses during genetic algorithm

Usage

## S3 method for class 'pso'
plot(x, y = NULL, ...)

Arguments

Value

This function does not return anything.

Print method for "gena" object

Description

Print method for "gena" object

Usage

## S3 method for class 'gena'
print(x, ...)

Arguments

Value

This function does not return anything.

Print method for "pso" object

Description

Print method for "pso" object

Usage

## S3 method for class 'pso'
print(x, ...)

Arguments

Value

This function does not return anything.

Summary for "gena" object

Description

Summary for "gena" object

Usage

## S3 method for class 'summary.gena'
print(x, ...)

Arguments

Value

This function returns x input argument.

Summary for "pso" object

Description

Summary for "pso" object

Usage

## S3 method for class 'summary.pso'
print(x, ...)

Arguments

Value

This function returns x input argument.

Particle Swarm Optimization

Description

This function allows to use particle swarm algorithm for numeric global optimization of real-valued functions.

Usage

pso(
  fn,
  gr = NULL,
  lower,
  upper,
  pop.n = 40,
  pop.initial = NULL,
  pop.method = "uniform",
  nh.method = "random",
  nh.par = 3,
  nh.adaptive = TRUE,
  velocity.method = "hypersphere",
  velocity.par = list(w = 1/(2 * log(2)), c1 = 0.5 + log(2), c2 = 0.5 + log(2)),
  hybrid.method = "rank",
  hybrid.par = 2,
  hybrid.prob = 0,
  hybrid.opt.par = NULL,
  hybrid.n = 1,
  constr.method = NULL,
  constr.par = NULL,
  random.order = TRUE,
  maxiter = 100,
  is.max = TRUE,
  info = TRUE,
  ...
)

Arguments

Details

Default arguments have been set in accordance with SPSO 2011 algorithm proposed by M. Clerc (2012).

To find information on particular methods available via pop.method, nh.method, velocity.method, hybrid.method and constr.method arguments please see 'Details' section of gena.population, pso.nh, pso.velocity, gena.hybrid and gena.constr correspondingly.

It is possible to provide manually implemented functions for population initialization, neighbourhoods creation, velocity updated, hybridization and constraints in a similar way as for gena.

By default function does not impose any constraints upon the parameters. If constr.method = "bounds" then lower and upper constraints will be imposed. Lower bounds should be strictly smaller then upper bounds.

Currently the only available termination condition is maxiter. We are going to provide some additional termination conditions during future updates.

Infinite values in lower and upper are substituted with -(.Machine$double.xmax * 0.9) and .Machine$double.xmax * 0.9 correspondingly.

By default if gr is provided then BFGS algorithm will be used inside optim during hybridization. Otherwise Nelder-Mead will be used. Manual values for optim arguments may be provided (as a list) through hybrid.opt.par argument.

For more information on particle swarm optimization please see M. Clerc (2012).

Value

This function returns an object of class pso that is a list containing the following elements:

• par - particle (solution) with the highest fitness (objective function) value.

• value - value of fn at par.

• population - matrix of particles (solutions) of the last iteration of the algorithm.

• counts - a two-element integer vector giving the number of calls to fn and gr respectively.

• is.max - identical to is.max input argument.

• fitness.history - vector which i-th element is fitness of the best particle in i-th iteration.

• iter - last iteration number.

References

M. Clerc (2012). Standard Particle Swarm Optimisation. HAL archieve.

Examples

## Consider Ackley function

fn <- function(par, a = 20, b = 0.2)
{
  val <- a * exp(-b * sqrt(0.5 * (par[1] ^ 2 + par[2] ^ 2))) +
         exp(0.5 * (cos(2 * pi * par[1]) + cos(2 * pi * par[2]))) -
         exp(1) - a
  return(val)
}

# Maximize this function using particle swarm algorithm

set.seed(123)
lower <- c(-5, -100)
upper <- c(100, 5)
opt <- pso(fn = fn, 
           lower = lower, upper = upper,
           a = 20, b = 0.2)
print(opt$par)


## Consider Bukin function number 6

fn <- function(x, a = 20, b = 0.2)
{
  val <- 100 * sqrt(abs(x[2] - 0.01 * x[1] ^ 2)) + 0.01 * abs(x[1] + 10)
  return(val)
}

# Minimize this function using initially provided
# position for one of the particles
set.seed(777)
lower <- c(-15, -3)
upper <- c(-5, 3)
opt <- pso(fn = fn, 
           pop.init = c(8, 2),
           lower = lower, upper = upper,
           is.max = FALSE)
print(opt$par)

Neighbourhood

Description

Constructs a neighbourhood of each particle using particular topology.

Usage

pso.nh(pop.n = 40, method = "ring", par = 3, iter = 1)

Arguments

Details

If method = "ring" then each particle will have par[1] neighbours. By default par[1] = 3. See section 3.2.1 of M. Clerc (2012) for additional details. If method = "wheel" then there is a single (randomly selected) particle which informs (and informed by) other particles while there is no direct communication between other particles. If method = "random" then each particle randomly informs other par[1] particles and itself. Note that duplicates are possible so sometimes each particle may inform less then par[1] particles. By default par[1] = 3. See section 3.2.2 of M. Clerc (2012) for more details. If method = "star" then all particles are fully informed by each other. If method = "random2" then each particle will be self-informed and informed by the j-th particle with probability par[1] (value between 0 and 1). By default par[1] = 0.1.

Value

This function returns a list which i-th element is a vector of particles' indexes which inform i-th particle i.e. neighbourhood of the i-th particle.

References

Maurice Clerc (2012). Standard Particle Swarm Optimisation. HAL archieve.

Examples

# Prepare random number generator
set.seed(123)

# Ring topology with 5 neighbours
pso.nh(pop.n = 10, method = "ring", par = 5)

# Wheel topology
pso.nh(pop.n = 10, method = "wheel")

# Star topology
pso.nh(pop.n = 10, method = "star")

# Random topology where each particle 
# randomly informs 3 other particles
pso.nh(pop.n = 10, method = "random", par = 3)

# Random2 topology wehere each particle could
# be informed by the other with probability 0.2
pso.nh(pop.n = 10, method = "random2", par = 0.2)

Velocity

Description

Calculates (updates) velocities of the particles.

Usage

pso.velocity(
  population,
  method = "hypersphere",
  par = list(w = 1/(2 * log(2)), c1 = 0.5 + log(2), c2 = 0.5 + log(2)),
  velocity,
  best.pn,
  best.nh,
  best.pn.fitness,
  best.nh.fitness,
  iter = 1
)

Arguments

Details

If method = "classic" then classical velocity formula is used:

v_{i,j,(t+1)}=w\times v_{i,j,t} + c_{1}\times u_{1,i,j} \times b^{pn}_{i,j,t} + c_{2}\times u_{2,i,j} \times b^{nh}_{i,j,t}

where v_{i, j, t} is a velocity of the i-th particle respect to the j-th component at time t. Random variables u_{1,i,j} and u_{2,i,j} are i.i.d. respect to all indexes and follow standard uniform distribution U(0, 1). Variable b^{pn}_{i,j,t} is j-th component of the best known particle's (personal) position up to time period t. Similarly b^{nh}_{i,j,t} is j-th component of the best of best known particle's position in a neighbourhood of the i-th particle. Hyperparameters w, c_{1} and c_{2} may be provided via par argument as a list with elements par$w, par$c1 and par$c2 correspondingly.

If method = "hypersphere" then rotation invariant formula from sections 3.4.2 and 3.4.3 of M. Clerc (2012) is used with arguments identical to the classical method. To simulate a random variate from the hypersphere function rhypersphere is used setting type = "non-uniform".

In accordance with M. Clerc (2012) default values are par$w = 1/(2 * log(2)), par$c1 = 0.5 + log(2) and par$c2 = 0.5 + log(2).

Value

This function returns a matrix which i-th row represents updated velocity of the i-th particle.

References

Maurice Clerc (2012). Standard Particle Swarm Optimisation. HAL archieve.

Hypersphere

Description

Simulates uniform random variates from the hypersphere.

Usage

rhypersphere(n, dim = 2, radius = 1, center = rep(0, dim), type = "boundary")

Arguments

Value

The function returns a vector of random variates.

Examples

set.seed(123)
# Get 5 random uniform variates from 3D hypersphere
# of radius 10 centered at (2, 3, 1)
rhypersphere(n = 5, dim = 3, radius = 10, center = c(2, 3, 1))

Summarizing gena Fits

Description

Summarizing gena Fits

Usage

## S3 method for class 'gena'
summary(object, ...)

Arguments

Value

This function returns the same list as gena function changing its class to "summary.gena".

Summarizing pso Fits

Description

Summarizing pso Fits

Usage

## S3 method for class 'pso'
summary(object, ...)

Arguments

Value

This function returns the same list as pso function changing its class to "summary.pso".