# Chaotic Honeybees Optimization Algorithms Approach for Traveling Salesperson Problem

Authors: Pedro Palominos; Carla Ortega; Miguel Alfaro; Guillermo Fuertes; Manuel Vargas; Mauricio Camargo; Victor Parada; Gustavo Gatica
Journal: Complexity
 (2022)

---

## Abstract

Due to the difficulty in solving combinatorial optimization problems, it is necessary to improve the performance of the algorithms by improving techniques to deal with complex optimizations. This research addresses the metaheuristics of marriage in honey-bees optimization (MBO) based on the behavior of bees. The current study proposes a technique for solving combinatorial optimization problems within proper computation times. The purpose of this study focuses on the travelling salesperson problem and the application of chaotic methods in important sections of the MBO metaheuristic. Three experiments were conducted to measure the efficiency and quality of the solutions: (1) MBO with chaos to generate initial solutions (MBO2); (2) MBO with chaos in the workers (MBO3); and (3) MBO with chaos to generate initial solutions and the workers (MBO4). The application of chaotic functions in MBO was significantly better at solving the travelling salesperson problem.


---

## Body

## 1. Introduction

The term metaheuristics was first introduced by Glover [1] in 1986. Establishing a mathematical model for complex systems is often a difficult task [2]. Classical methods depend on the type of objective function and constraint and on the type of variable used in modeling the problem. In addition, the effectiveness of classical algorithms is highly dependent on the solution section (convex or nonconvex), the number of decision variables, and the number of constraints in the modeling of the problem [3]. However, classical optimization algorithms are insufficient for large-scale combinational and nonlinear problems [4]. According to Bingol and Alatas [5], artificial intelligence-based metaheuristic methods are generally categorized as physics-based, biologically-based, social-based, musical-based, sport-based, swarm-based, plant-based, chemistry-based, light-based, mathematics-based, and water-based. Efficient hybrid methods have also been proposed by combining these algorithms. Based on a single solution (direct search algorithms), the following algorithms can be found: simulate annealing (SA) [6], taboo search (TS) [7, 8], random walk (RW) [9], and hill climbing (HC) [10], among others, and population-based algorithms such as spider monkey optimization (SMO) [11]; particle swarm optimization (PSO) [12–14]; ant colony optimization (ACO) [15]; artificial immune system (AIS) [16]; whale optimization [17]; genetic algorithm (GA) [18, 19]; firefly algorithm [20]; grey wolf optimizer (GWO) [21]; bee algorithm (BA) [22]; artificial bee colony (ABC) [23]; queen bee evolution (QBE) [24]; bee system (BS) [25, 26]; bee colonies optimization (BCO) [27]; BeeAdHoc [28, 29]; and marriage in honey bees optimization (MBO) [30, 31] and its different versions such as honey bees mating optimization (HBMO) [32, 33] fast marriage in honey bees optimization (FMHBO) [34], and honey bees optimization (HBO) [35].MBO is a metaheuristic proposed by Abbass [30] in 2001, where queens, drones, larvae, and workers interact to solve different problems. The genetic material of queens and drones gives rise to larvae cared for by the workers, and the latter heuristics improve the solution [36]. According to Shamsaldin et al. [37], MBO has advantages over GA when performing a local search by iteration.Alfaro et al. [38, 39] researched dynamic systems to find new heuristic techniques that generated better solutions or improved the convergence speed. They proposed using chaos theory as an optimization tool. Chaos is a nonlinear deterministic phenomenon, sensitive to initial conditions and capable of determining all states within a range without repetition [40, 41]. It is used in optimization because it is easy to traverse the search space and avoids local minimums.Other studies [42–46] related to evolutionary algorithms enhanced by chaotic maps have been described. In those algorithms, the generation of random values for the different parameters in the model is replaced by the application of chaotic maps for generating those parameters. The performance of these metaheuristics with evolutionary operators generated with chaotic numbers is better than that with evolutionary operators with random numbers. Evolutionary algorithms reduce premature convergence and increase performance optimization. Numerous chaotic systems, such as the Henon map, Logistic map, Rossler map, Tent map, and Zaslavskii map, can be used. The Lorenz system results in the most efficient one [47].Fuertes et al. [47] have studied the performance of chaotic numbers used with a raw genetic algorithm. Their work proposes a new edge for the chaotic genetic algorithm (CGA) and the importance of entropy in the initial population. CGA uses chaotic maps to modify the stochastic parameters of the genetic algorithm. The algorithm modifies the parameters of the initial population using chaotic series and then analyzes the entropy of such a population. The numerical experiment demonstrates a correlation between entropy and the performance of the algorithm. The study concludes that the chaotic series with larger entropy used in CGA has better performance optimization than random series.Vargas et al. [48] studied the effect of entropy on the performance of the genetic algorithm modified using a time series of earthquakes and winds. In this study, the genetic algorithm modifies the stochastic parameters of the genetic algorithm with the chaotic sequences of numbers from the earthquake and wind time series. The experiment demonstrates better optimization performance compared to the stochastic genetic algorithm.The purpose of this study is to analyze and evaluate whether chaos theory improves solution quality in MBO for the travelling salesperson problem (TSP) using chaotic methods in essential MBO sections. This research measured this technique’s efficiency and solution quality.To measure the obtained result, all experiments will be carried out using TSP. TSP is difficult to solve and requires a lot of computational time to find optimal solutions, so most prefer to use heuristic techniques [49–51]. Therefore, TSP was chosen to measure the efficiency of chaos in MBO.


### 1.1. Related Literature

Pan et al. [52] proposed a new neuronal circuit with a synapse module, a neural module, and chaotic characteristics. The study designed and studied a chaotic neural network (CNN), applied associative memory to bipolar images, and verified that the circuit was capable of chaotic associative memory.Lin [53] proposed a predefined time stability theorem based on the Lyapunov function to facilitate the antisynchronization of the two CNNs. The engineered driver causes antisync errors to converge within a predefined period.Long et al. [54] studied finite-time control theory for a class of state-based switched inertial chaotic neural networks (SBSCINN) with distributed delays in practical applications. These require fast or limited convergence time to ensure CNN stability. The study also generated new criteria for stabilizing SBSCINNs.Chen et al. [55] studied a fractional order discrete Hopfield neural network (FODHNN) in an image encryption system, showed that it had complex nonlinear dynamic behavior, and designed a timing controller.Srivastava et al. [56] proposed an algorithm to improve data security using a hybrid model that uses an adaptive coding technique with a CNN hopfield. Their results indicated that the speed and security of the transmitted data are better than those of traditional algorithms. Additionally, the calculation time for the proposed algorithm was shorter than that of traditional algorithms.Jiang and Chen [57] proposed a blind restoration model for super-resolution images based on CNN. The average values of blind restore time, image sharpness, and power consumption were 9.273 ms, 99.045%, and 118.524 J, respectively. The model restoration performance was good, and the super-resolution blind image was the most like the original image.Lian et al. [58] proposed a new discontinuous Lyapunov function to solve local synchronization control for CNN with sampled data and saturation actuators. Additionally, the authors provided a method for optimizing the design of sampled-data controllers. Their proposed method could be extended to discrete-time CNN by using finite-sum inequalities.Lee [59] proposed the use of CNN in quantum finance. They studied quantum finance oscillators (QFOs) using Lee oscillators and different QFO applications, including quantum financial prediction using CNN, deep learning CNN, and chaotic multi-agent based smart trading systems. They reported that chaotic neural oscillators could be adopted in almost any artificial intelligence tool or technology to improve efficiency and performance.Lahmiri and Bekiros [60] implemented deep learning CNN to forecast the price of the three most traded digital currencies, Bitcoin, Digital Cash, and Ripple. The time series of all digital currencies exhibits fractal dynamics, long memory, and self-similarity. They found that deep learning forecasted the inherent chaotic dynamics of cryptocurrency markets very efficiently.Wang and Zhang [61] designed two fractional observers for CNN with and without parametric uncertainties. Their observation errors eventually converged in a small region using fractional stability criteria, and their proposed observers were robust.Chen et al. [62] proposed a new multiheating dynamic PSO framework that combines the chaotic grouping mechanism (CGM) and the dynamic regrouping strategy (DRS), which they named CGPSO-DRS. Their method improved the convergence speed and time to find the global optimum compared to similar population-based approaches and next-generation PSO variants.Rosalie et al. [63] proposed a Bayesian optimization model to select Rössler system parameters used in chaotic ant colony optimization for coverage (CACOC). They developed this model to manage the movements of an unmanned aerial vehicle swarm (UAVs), and their approach improved the UAV swarm coverage speed.Tummala et al. [64] proposed a new algorithm called “war strategy optimization” (WSO) where the war strategy is modeled as an optimization process where each soldier moves dynamically toward the optimal value. The performance of the algorithm is compared with the other ten popular metaheuristic algorithms. Experimental results for several optimization problems demonstrate the superiority of WSO.Alatas [65] developed a computational method inspired by the types and occurrence of chemical reactions called the artificial chemical reaction optimization algorithm (ACROA). The method is based on the principles of chemical reactions. These chemical reactions transform reactants through a sequence of reactions into products. The initial reactants are distributed in the feasible search region, obtaining high-quality partial solutions and achieving optimal solutions in a few cycles. This method includes global and local search capabilities and does not require a local search method to refine the search.The proposed Chaotic MBO technique is based on a well-known metaheuristic-type method with more than a decade of research. The calibration of the MBO metaheuristic has been extensively studied, and it allows us to start the research from a very robust state-of-the-art. The research took advantage of this knowledge to implement chaotic determinism in the metaheuristic and to search for improvements in the optimization technique.Compared to other techniques presented by many authors, the chaotic deterministic characteristics of the phenomena used for this study, such as chemical processes, human behavior, chaotic maps, or quantum phenomena, establish a nonlinear search for viable solutions. The results of these techniques have presented sustainably robust and high-quality solutions. Therefore, the selection of the MBO metaheuristic and its modification with chaotic determinism allows us to quickly obtain results of high accuracy.


### 1.2. Limitations and Contributions

The development of this new computational method was inspired by the behavior of honeybees through the interaction of the different members of the hive. This method establishes the global search capability through crossover cycles of the best partial solutions. The proposed research improves the performance of the metaheuristic methodology by applying chaotic determinism for creating drones and the process of worker bees. The workers improve the solution quality with chaotic modifications like chaos simulated annealing (CSA), chaotic tabu search neural network based on a path (CTSNN), or chaotic swap local search (CSLS). Previous research has demonstrated the robustness of the MBO method for discarding partial optimal solutions; the proposed modification improves the ability to search for optimal solutions.The metaheuristic method was used to solve the TSP problem; the results show an improvement in the performance of the chaotic MBO metaheuristics for this problem. For other engineering problems, it is necessary to replicate the experiment with new parameters.Finally, this metaheuristic method does not require many parameters for its implementation; the search time does not increase quite significantly; and it can be adapted for multiobjective optimization models or parallel distributions.


## 1.1. Related Literature

Pan et al. [52] proposed a new neuronal circuit with a synapse module, a neural module, and chaotic characteristics. The study designed and studied a chaotic neural network (CNN), applied associative memory to bipolar images, and verified that the circuit was capable of chaotic associative memory.Lin [53] proposed a predefined time stability theorem based on the Lyapunov function to facilitate the antisynchronization of the two CNNs. The engineered driver causes antisync errors to converge within a predefined period.Long et al. [54] studied finite-time control theory for a class of state-based switched inertial chaotic neural networks (SBSCINN) with distributed delays in practical applications. These require fast or limited convergence time to ensure CNN stability. The study also generated new criteria for stabilizing SBSCINNs.Chen et al. [55] studied a fractional order discrete Hopfield neural network (FODHNN) in an image encryption system, showed that it had complex nonlinear dynamic behavior, and designed a timing controller.Srivastava et al. [56] proposed an algorithm to improve data security using a hybrid model that uses an adaptive coding technique with a CNN hopfield. Their results indicated that the speed and security of the transmitted data are better than those of traditional algorithms. Additionally, the calculation time for the proposed algorithm was shorter than that of traditional algorithms.Jiang and Chen [57] proposed a blind restoration model for super-resolution images based on CNN. The average values of blind restore time, image sharpness, and power consumption were 9.273 ms, 99.045%, and 118.524 J, respectively. The model restoration performance was good, and the super-resolution blind image was the most like the original image.Lian et al. [58] proposed a new discontinuous Lyapunov function to solve local synchronization control for CNN with sampled data and saturation actuators. Additionally, the authors provided a method for optimizing the design of sampled-data controllers. Their proposed method could be extended to discrete-time CNN by using finite-sum inequalities.Lee [59] proposed the use of CNN in quantum finance. They studied quantum finance oscillators (QFOs) using Lee oscillators and different QFO applications, including quantum financial prediction using CNN, deep learning CNN, and chaotic multi-agent based smart trading systems. They reported that chaotic neural oscillators could be adopted in almost any artificial intelligence tool or technology to improve efficiency and performance.Lahmiri and Bekiros [60] implemented deep learning CNN to forecast the price of the three most traded digital currencies, Bitcoin, Digital Cash, and Ripple. The time series of all digital currencies exhibits fractal dynamics, long memory, and self-similarity. They found that deep learning forecasted the inherent chaotic dynamics of cryptocurrency markets very efficiently.Wang and Zhang [61] designed two fractional observers for CNN with and without parametric uncertainties. Their observation errors eventually converged in a small region using fractional stability criteria, and their proposed observers were robust.Chen et al. [62] proposed a new multiheating dynamic PSO framework that combines the chaotic grouping mechanism (CGM) and the dynamic regrouping strategy (DRS), which they named CGPSO-DRS. Their method improved the convergence speed and time to find the global optimum compared to similar population-based approaches and next-generation PSO variants.Rosalie et al. [63] proposed a Bayesian optimization model to select Rössler system parameters used in chaotic ant colony optimization for coverage (CACOC). They developed this model to manage the movements of an unmanned aerial vehicle swarm (UAVs), and their approach improved the UAV swarm coverage speed.Tummala et al. [64] proposed a new algorithm called “war strategy optimization” (WSO) where the war strategy is modeled as an optimization process where each soldier moves dynamically toward the optimal value. The performance of the algorithm is compared with the other ten popular metaheuristic algorithms. Experimental results for several optimization problems demonstrate the superiority of WSO.Alatas [65] developed a computational method inspired by the types and occurrence of chemical reactions called the artificial chemical reaction optimization algorithm (ACROA). The method is based on the principles of chemical reactions. These chemical reactions transform reactants through a sequence of reactions into products. The initial reactants are distributed in the feasible search region, obtaining high-quality partial solutions and achieving optimal solutions in a few cycles. This method includes global and local search capabilities and does not require a local search method to refine the search.The proposed Chaotic MBO technique is based on a well-known metaheuristic-type method with more than a decade of research. The calibration of the MBO metaheuristic has been extensively studied, and it allows us to start the research from a very robust state-of-the-art. The research took advantage of this knowledge to implement chaotic determinism in the metaheuristic and to search for improvements in the optimization technique.Compared to other techniques presented by many authors, the chaotic deterministic characteristics of the phenomena used for this study, such as chemical processes, human behavior, chaotic maps, or quantum phenomena, establish a nonlinear search for viable solutions. The results of these techniques have presented sustainably robust and high-quality solutions. Therefore, the selection of the MBO metaheuristic and its modification with chaotic determinism allows us to quickly obtain results of high accuracy.


## 1.2. Limitations and Contributions

The development of this new computational method was inspired by the behavior of honeybees through the interaction of the different members of the hive. This method establishes the global search capability through crossover cycles of the best partial solutions. The proposed research improves the performance of the metaheuristic methodology by applying chaotic determinism for creating drones and the process of worker bees. The workers improve the solution quality with chaotic modifications like chaos simulated annealing (CSA), chaotic tabu search neural network based on a path (CTSNN), or chaotic swap local search (CSLS). Previous research has demonstrated the robustness of the MBO method for discarding partial optimal solutions; the proposed modification improves the ability to search for optimal solutions.The metaheuristic method was used to solve the TSP problem; the results show an improvement in the performance of the chaotic MBO metaheuristics for this problem. For other engineering problems, it is necessary to replicate the experiment with new parameters.Finally, this metaheuristic method does not require many parameters for its implementation; the search time does not increase quite significantly; and it can be adapted for multiobjective optimization models or parallel distributions.


## 2. Research Methodology

### 2.1. The Travelling Salesman Problem

TSP is a well-studied problem in combinatorial optimization because it is difficult to solve. A vendor wants to visitn cities only once and finally return to the city where they started their journey. The idea is to minimize the journey cost, where the cost of a city i to a city j corresponds to the distance travelled between said cities (Cij). Computer scientists call it an NP-complex problem [66], because it cannot be solved by an algorithm in polynomial time.Mathematically, TSP can be expressed as follows:Objective function:(1)min∑i,jcij⋅xij.Subject to:(2)∑ixij=1∀j≠i,∑jxij=1∀i≠j,u1=1,(3)2≤ui≤n∀i≠1,(4)ui−uj+1≤n−11−xij∀i≠1,∀j≠1,(5)ui≥0∀i,xij∈0,1∀i,j.Constraints (2)–(5) reach together to eliminate the subpaths in the solution.TSP can be used to solve different problems. Robotics can minimize the number of machine movements, where each movement has an associated cost [67]. Distribution logistics can find the minimum path from the point of origin, passing through known places and finally returning to the point of origin [68]. This happens in a dispatch system that minimizes fuel usage and ensures that the products are delivered as quickly as possible.


### 2.2. Honeybee Optimization Metaheuristics

Among the methods inspired by nature is the artificial bee colony (ABC) algorithm, which belongs to the swarm intelligence family. Within the bee colony methods, there are two distinct types of metaheuristics. The first approach is the bee collecting pollen algorithm (BCPA), where bees seek to obtain the maximum amount of nectar [69].MBO is based on the social behavior of bees in which queens, drones, larvae, and workers interact to solve overly complex combinatorial problems.


#### 2.2.1. Description of the Metaheuristics

A bee colony can be formed in two different ways. The first method is called independent foundation, in which the colony is created by one or more queens that build a nest, produce eggs, and feed the larvae. Later, the queens only take charge of laying eggs, and the workers conduct the work in the hive. The second method is called swarming, in which one or more queens and workers create the colony from the original colony, and the division of labor begins immediately; therefore, the queens are only responsible for laying eggs. Colony formation by a single queen is called haplometrosis; otherwise, it is called pleometrosis. Additionally, if the hive has only one queen during its life cycle, it is monogynous, and if it has many queens, it is polygynous [70].Colonies are made up of queens, drones, larvae, and workers. Queens lay eggs, and male drones are haploid and can only affect the queen’s genotype via mutations. Fertilization occurs in the air, and the queen interbreeds 7 to 20 drones in one session. The drone’s genetic material is stored in the queen’s spermatheca, and a random mixture of the stored sperm fertilizes the queen’s eggs. The workers who care for the larvae also occasionally lay eggs.Our model randomly creates an initial population of queens, and its solutions are improved using workers. Each queen initiates a mating flight with initial energy and speed.The mating flight is a set of transitions where the queen moves in the air and mates (or does not mate) with a drone according to a probabilistic function. If mating occurs, the drone’s sperm is stored within the queen’s spermatheca. The flight ends when the queen’s energy has been exhausted or the spermatheca is at maximum capacity. After that, her energy and speed decreased.When the queen is in the nest, the drone’s genetic material is selected to create larvae cared for by the workers. The workers are represented by heuristics that seek to improve the larva (solution).Eventually, the best larvae become queens and replace the worst queens. The remaining queens and larvae are eliminated to begin a new flight cycle.The probability that queens and drones cross is represented by the annealing function, which depends on the queen’s fitness, the drone’s fitness, and the queen’s speed which is given as follows:(6)pq,d=mín1,e−lq,d/sq,where l (q, d) = dist(f(q), f(d)) corresponds to the absolute difference between the fitness of the queen and the drone, s(q) is the speed of the queen at that instant.The probability of mating is greater at the beginning of the flight since the queen’s energy is greater, and it is more likely the drone’s fitness is equivalent to the queen’s.As mentioned above, after each transition, speed and energy must be reduced. The speed must be adjusted according to the equation:(7)st+1=α⋅st,where α is a factor in the interval [0, 1].When energy is reduced, the following equation is used:(8)et+1=et−step,where step corresponds to the energy reduction, whose value is calculated according to the following equation:(9)step=0,5⋅e0M.M is the capacity of the spermatheca.Figure1 shows the five main stages of MBO: (1) mating flight in which the queen probabilistically selects the drones she will mate with and adds their sperm to her spermatheca; (2) creation of larvae by crossing between the genes of the queen and the drones; (3) using workers (heuristics) to improve the solution represented by the larvae; (4) update the worker’s fitness according to how well they improve the larvae; and (5) replace the worst queens with the best larvae.Figure 1
Genotype example.Abbass [31, 71] only considered a single queen and a single worker to solve a particular class of the satisfaction problem (3-SAT) and compared the results obtained with the WalkSAT and GSAT heuristics. MBO with the worker WalkSAT was more satisfying than WalkSAT alone, and so was MBO with the GSAT worker. The whole experiment was developed with parameters taken from biology.Teo and Abbass [72] changed the accepted trajectory between the drone and the queen. This time, the drone was accepted, and its sperm was added to the queen’s spermatheca if it was fitter than the queen or fulfilled the new annealing function which is given as follows:(10)pq,d=eΔf/sq,where Δ(f) corresponds to the difference between the fitness of the new trajectory and the queen and s(q) is the speed of the queen at that instant.Later, Abbass and Teo [36] changed the annealing function. Instead of using the difference between the queen’s fitness and the drone, they used the difference between the drone trajectories’ fitness to ensure that a trajectory was only accepted if it was better than the previous ones used. The new annealing function was determined by the following equation:(11)pq,d=e−ft−ft−1/sq,where (f(t) − f(t − 1)) corresponds to the difference in the fitness of the current trajectory and the previous trajectory and s(q) is the speed of the queen at that instant.The study addressed the 3-SAT problem using five workers: GSAT, random walk, probabilistic greedy, one-point crossover, and WalkSAT. The tests were conducted by taking each heuristic separately, using MBO with each worker, and using the workers together in MBO. The results showed that MBO behaves better when all the workers are used together and that the new modification [36] provides better results than the original version [31, 71] and the one proposed by [72].


#### 2.2.2. State of the Art MBOs

Phu-Ang and Thammano [73] proposed a new MBO-based memetic algorithm to solve the flexible workshop scheduling problem. They introduced four new features to the standard MBO algorithm to move the search away from the local optimum. Their proposed algorithm maintained two populations: the female and male, and their performance was compared with the dispersion search algorithms using path linkage, hybrid differential evolution, hybrid harmony search algorithms, and hybrid GA. MBO outperformed others on large and complex problem cases.Wen et al. [74] proposed a modified HBMO (MHBMO) to integrate process planning and scheduling (IPPS) with uncertain processing time and due date based on the fuzzy set by designing an uncertainty measurement objective calculation method. Their method solves the IPPS of multiple objectives effectively.Palominos et al. [75] proposed a new MBO metaheuristic to deal with shop flow scheduling problems by designing different mating flight spaces and achieved excellent results over the 120 instances they tested.Solgi and Loáiciga [76] identified seven algorithms for solving continuous optimization problems and showed that ABC, bee evolution for genetic algorithms (BEGA), and MBO were the most efficient.Vakil-Baghmisheh and Salim [77] proposed a modified version of MFMHBO and compared it with ABC, QBE, and FMHBO in four reference functions for various variables up to 100. MFMHBO was faster than the others in most cases, especially for the Griewank and Schwefel functions, and increased precision and the number of variables. In general, their modified algorithm is comparable to the other algorithms.Prabhakar and Lee [78] studied different MBO techniques for transformation-based three-level character selection using wavelets for prostate cancer classification.Çelik and Ülker [79] tried to solve TSP using MBO. TSP had previously been solved using heuristic methods (GA and SA), but their MBO algorithm converged toward the optimal solution in fewer iterations.Nayak and Naik [80] presented a hybrid metaheuristic pi-sigma neural network (PSNN) based on honey bee mating that successfully solved the data mining classification problem. Their approach combined HBMO with PSNN and compared it with other techniques such as GA, differential evolution (DE), and PSO. Their approach was stable and reliable and provided better classification precision than others.Yin et al. [81] designed an accurate and stable MBO algorithm for programming a scorching machine to minimize the weighted sum of work completion times.Zanbouri and Jafari Navimipour [82] proposed creating cloud computing services where services were selected based on trust and quality (QoS). They used the media grouping algorithm k to reduce the search space and HBMO to achieve global optimization. Their method worked efficiently in terms of computational time, producing more reliable services, although it was not efficient for large-scale data sets.Arun and Vijay Kumar [83] proposed a modified MBO for a view selection algorithm (MBOVSA) to select views with the lowest total cost and information processing time in a large data warehouse. MBOVSA selected better quality views than HRUA, one of the important greed-based view selection algorithms.


### 2.3. Experimental Settings

#### 2.3.1. TSP Modelling using MBO

Our technique is based on [30]. We implemented Abbass’ proposal with some modifications in TSP. The MBO metaheuristic is described in Algorithm 1.Algorithm 1: MBO Metaheuristics.
Define the number of queensQ, Workers W, larvae B, mating flights G, and spermatheca sizeDefine the energy and speed of the queensAssign each worker with a different heuristicAssign the genotype of the queens randomlyChoose a worker randomly and apply its genotype to the queen’s genotypeWhile the number of mating flights<= GFor queen = 1 throughQAssign queen power and speedAssign step = 0.5 × energy/MAssign path randomlyAs long as energy> 0 and the spermatheca is not fullThe queen moves between different states and probabilistically chooses dronesIf a drone is selected, thenAdd your sperm to the queen’s spermathecaEndUpgrade energy and speed of the queenWhileEnd ForGenerate larva by crossing and mutationUse workers to improve larvaUpdate the fitness of the workersWhile the best of larvae is better than the worst of queensReplace the worst queen with the best larvaRemove the best larva from the larvae listEnd whileRemove all larvaEnd While


#### 2.3.2. Defining Parameters

The proposed algorithm requires a series of parameters: (1) number of queen bees; (2) number of drones; (3) capacity of the queen’s spermatheca; (4) number of flights; (5) number of larvae; (6) initial energy; and (7) initial speed of the queen. Particular care was taken when choosing (7) (defined experimentally) because the speed influences the annealing function, which selects the genetic material of the drone. When (7) is low, it is more likely that a drone will not be chosen for a flight, and if it is high, it is more likely that all drones will be accepted. Additionally, since the fitness range can be broad in TSP, one should not select a fixed speed for all instances to avoid the above-mentioned problems. Therefore, (7) must be calculated considering the difference between the best and the worst queen when starting the algorithm. For parameter (6), values from previous studies published in [30] were used for the energy reduction and speed factors. Table 1 shows the parameter values used for the proposed MBO algorithm.Table 1
Algorithm parameters.
ParametersSource of obtainingValueNumber of queen beesExperimental3Drone quantityExperimental200Spermatheque capacityExperimental50Number of flightsExperimental100Number of larvaeExperimental30Initial energy[84]0.9Initial velocityExperimentalDifference between the best and the worst queen at the start of the algorithm.Energy reduction factor[30]step=0.5⋅e0/MSpeed reduction factor[30]α = 0.9


#### 2.3.3. Representing the Solution

The genotype stereotype will be the same for queens, drones, and larvae and corresponds to a sequence of cities to be visited by the traveler. The genotype length is the number of cities that must be stopped by and varies according to the problem. Each gene contains the cities and the corresponding order. For example, if there is a three in position four, the city three will be after the fourth.Figure2 shows an example solution for a problem with six cities where the optimal order must be 1, 3, 6, 2, 5, and 4.Figure 2
HBMO algorithm [71].Each genotype has an associated fitness that corresponds to the tour distance.


#### 2.3.4. Generating an Initial Population of Bees, Queens, and Drones

The initial generation of queens can be selected in three diverse ways. The first two queens are generated on the basis of two tour construction heuristics for TSP: the nearest neighbor [85] and the closest insertion [86]. The remaining queens and drones are randomly generated to introduce diversity [36].


#### 2.3.5. Selection of Drone Genetic Material

Drone genetic material is selected by the MBO annealing function [30]. Drone (d) selection is more likely when its fitness is like the queen (q) and when the queen begins her flight.A random number between 0 and 1 is generated to select the drone. When this number is less than or equal to the one provided by the annealing function, the drone’s genetic material is added to the queen’s spermatheca (if the random number ≤pq,d, then add genetic material from the drone to the queen’s spermatheca).


#### 2.3.6. Mating

After completing a flight, each queen returns to the nest to lay eggs. When all the queens have returned, the genetic materials used to generate larvae are selected. The number of larvae to be generated is an algorithm parameter, and each larva can come from any queen and any drone possessed by the chosen queen.Queens are chosen by fitness, i.e., fitter queens are more likely to be selected. Drone genetic material is picked up randomly, and larvae are generated by crossing.The crossing steps are as outlined below:(1)
Drones can only contribute half their genetic material because they are haploid [30, 87, 88]. Half of the drone’s genes tag like this to ensure.(2)
A random interval and segment (i, j), where i < j. i represents the beginning of the segment, and j represents the end.(3)
Each queen’s genes are checked, and they will be selected if they were not considered in the marked drone’s interval. The set of genes chosen by the queen generates a list.(4)
The first selectedi − 1 queen genes become part of the first i − 1 genes of the larva.(5)
The unlabeled genes in the marked drone segment (i, j) become part of the larva at the same position as the original gene.(6)
The queen’s gene that has not been considered goes back in the same position as the larva.(7)
From positionj + 1 the larva contains all other non-considered queen genes in the same order.Figure3 shows an example of the crossing process between a drone and a queen.Figure 3
Crossing.


#### 2.3.7. Larva Mutation

The larva mutates after being formed by the following steps:(1)
Each gene is visited, and a random number between 0 and 1 is generated.(2)
The gene is marked if the random number is less than or equal to a mutation probability.(3)
Marked genes are mutated.We defined the mutation probability defined as 1% per [71].


#### 2.3.8. Breeding Larvae Using Workers

We considered three types of workers corresponding to SA [89], path-based tabu search neural network (TSNN) [90], and swap local search (SLS) [91]. TSP created crossovers in our solutions, so we included an additional worker called the 2-opt [92] heuristic to avoid them. The assignment of this worker offers the solutions provided by the other three workers who were not chosen.


#### 2.3.9. Updating the Queen Bee Population

The queen population is updated after a flight is finished. We use the same process as the original MBO algorithm (Algorithm2).Algorithm 2: Update the original queen population.
While the best of larvae is better than the worst of queensReplace the worst queen with the best larvaRemove the best larva from the larvae listEnd whileThe solutions with the fittest queens and larvae generate new queens. The remaining solutions are discarded but can be considered in the next iteration. A new flight begins after the queen population has been updated, unless the stop criterion is met.


#### 2.3.10. Chaotic MBO Modeling

Chaos can be a useful tool for finding better solutions, but which MBO processes would benefit from chaos theory?The parts of the MBO metaheuristic (Algorithm1) must be considered, and a proposal to introduce chaos must be defined to answer this question. Chaos has successfully created initial solutions in other heuristics [93, 94].Drones are generated like queens, for which chaos better explores the solution space.In MBO, the workers improve the solutions. The effect of chaos cannot be evaluated if new chaotic heuristics replace MBO ones. Instead, chaos must be introduced into specific parts of the workers.We carried out different types of tests to measure the performance of chaos in the different MBO stages: (1) tests with traditional MBO; (2) tests with the use of chaos only in the generation of initial solutions (queens and drones); (3) tests with the use of chaos only in the workers; and (4) tests with the use of chaos in the generation of initial solutions and the workers.


#### 2.3.11. Chaotic Generation of Initial Queen and Drone Populations

The initial populations of queen bees and the generation of the drone population are generated using chaos. One queen is born using the nearest neighbor heuristic, and another uses the closest insertion heuristic. There are two heuristics in both versions to ensure good solutions and improve the algorithm’s performance.The remaining queens are born using chaotic maps from a randomly initialized solution vector where each gene takes a value between 0 and 1 for each value within the logistic map solution as shown in the following equation:(12)Xn+1=rXn1−Xn,where r is the control parameter (in this case, r = 4), whose value ensures that the system will enter a state of chaos. r has values between 0 and 1, except for the fixed points of 0.25, 0.5, and 0.75.When a queen starts her flight, the first drone drops randomly and has a value between 0 and 1. The remaining drones are on the logistic map the same way as the queens.


#### 2.3.12. Chaotic Breeding of Larvae using Workers

After generating the queens and larvae by crossing and mutations, the workers improve the solution quality.We use the same three workers with chaotic modifications (chaos simulated annealing (CSA), chaotic tabu search neural network based on a path (CTSNN), and chaotic swap local search (CSLS). Each worker performs a job that leads to better solutions.We included a chaotic version of the traditional MBO to avoid crossing TSPs in the solutions. This 2-opt heuristic uses the output of any chaotic worker mentioned above. The chaotic workers use the same fitness value at the beginning of the algorithm.


#### 2.3.13. Chaos Simulated Annealing (CSA) Heuristic

The first chaotic worker follows the same idea of simulated annealing (SA) in the traditional MBO modeling but with chaotic alterations (CSA).CSA was created by Mingjun and Huan Wen [95] to take advantage of SA in combinatorial optimization problems and the improvements it brings to chaos theory.The problem:(13)Pminfx,xi∈ai,bi,i=1,2,…,n.CSA and SA have the same structure but differ in two ways. The first difference is when creating an initial solutionx0 = (x01, x02,…, x0n) using the logistic map to obtain zki, according to following equation:(14)x0i=ai+bi−ai⋅zki,i=1,…,n.The second difference is when selecting a neighboring solutionym = (xm1, xm2,…, xmn), where one of the current solution variables is random and it´s modified using the following logistic map again:(15)ymi=xmi+α⋅bi−ai⋅≤zkm,where α is a chaos weighting parameter.Algorithm3 shows the CSA heuristic proposed by [95].Algorithm 3: Pseudocode for CSA.
Chaotic initialization, giving an initial value toz0Generate the chaotic variableszki, i = 1, 2,…, n, using the logistic functionThe initial solutionx0 = (x01, x02,…,x0n) is produced by the formula x0i = ai + (bi − ai) × zki, i = 1,...,n.Initialize maximum temperatureTmax, minimum temperature Tmin and the number of iterations LmaxT =  Tmaxm = 0Set as the best solutionx∗=x0Set as the best fitnessf∗=fx∗While (T > Tmin)While (m ≤ Lmax)Generate a new solutionym = (xm1, xm2,…,xmn) based on the formula ymi = xmi + α × (bi − ai) × zkm, picking out i randomly and zkm is delivered by the logistics mapα = α × e-βΔE∗  = f(ym) − f∗ΔE = f(ym) − f(xm)If (ΔE∗ ≤ 0)x∗  =  ymf∗=fx∗If (ΔE ≤ 0)xm+1 =  ymf∗=fx∗If (ΔE > 0)Accept status with probabilitye−ΔE/Tm = m + 1End whileLmax  =  Lmax + d, m = 0T = δ × TEnd whileDeliver the best solution foundx∗ and best fitness f∗We worked on chaotic genotypes where each gene has a chaotic variable between 0 and 1. However, not all solutions have a chaotic genotype, such as the larvae born after applying the crossover and mutation operators, which generate a solution based on the traditional MBO. The genotype is a chaotic solution.Furthermore, chaotic MBO uses the same parameters as traditional MBO. Still, it is necessary to define the value of the chaos weighting (α), which we determined as one by experimentation.It is also essential to define the chaotic function used, which is determined by the original CSA study and corresponds to the logistic map:(16)zk+1=4zk1−zk.


#### 2.3.14. Path-Based Chaotic Tabu Search Neural Network Heuristic

The second worker used for traditional MBO is path-based TSNN. Aono et al. [96] proposed a method based on the studies of [90], which incorporates chaotic noise into the equations.The new algorithm also relies on 2-opt moves to generate its TSP solutions. The new chaotic noise procedure is based on equations (17)–(19).(17)ξijt+1=βΔijt,(18)ζijt+1=−α∑d=0s−1krdxijt−d,(19)xijt+1=ξijt+1+ζaiajt+1+ζijt+1+γzijt+1,where Δij is the difference between the fitness of the current solution (D0) and the solution with a 2-opt movement between cities i and j (Dij): Δij (t) =  D0(t)  − Dij(t), β is the weight given to Δij(t), α is the weighting parameter of the taboo effect, kr is the decay parameter of the taboo effect, xij(t) is the output of neuron i, j in time t, and γ is the weight given to the chaotic noise zij.According to Aono et al. [96], the logistic map should generate zij in a chaotic form, as shown in the following equation:(20)zijt+1=3.828zijt1−zijt.The same parameters can be used for both the path-based TSNN and its chaotic version (CTSNN). The tests show different values ofγ in the interval [20, 30000], generating the best results for γ = 10000.


#### 2.3.15. Chaotic Swap Local Search Heuristic

Like the chaotic heuristic described above, chaotic noise gets the worker to swap local search and creates the chaotic swap local search.To develop this new method, all swap movements (i, j), where (i, j) = (j, i) e i ≠ j, as in the following equation:(21)ξijt+1=Δijt+γzijt,where Δij is the difference between the values of the current solution’s objective function (D0) and the solution with a swap movement of cities i and j (Dij): Δij (t) =  D0(t)  − Dij(t) and γ is the weight given to the chaotic noise zij.The next movement selected corresponds to the one with the highest value ofξij, which is not necessarily the one that delivers the best fitness because of the chaotic variable, but the one with the swap movement between cities i and j. Finally, the solution with the best fitness is generated when all iterations are carried out.This worker uses the same parameters as the nonchaotic version, adding chaotic noise through the following equation:(22)zijt+1=3.828zijt1−zijt.We selected the chaos weighting experimentally, taking values ofγ in the interval [20, 100] and obtaining the best results for γ = 40.


#### 2.3.16. Chaotic MBO Flow Chart

Figure4 shows a flow chart describing the chaotic MBO metaheuristic procedure described above to solve TPS.Figure 4
Chaotic MBO algorithm flow chart.


## 2.1. The Travelling Salesman Problem

TSP is a well-studied problem in combinatorial optimization because it is difficult to solve. A vendor wants to visitn cities only once and finally return to the city where they started their journey. The idea is to minimize the journey cost, where the cost of a city i to a city j corresponds to the distance travelled between said cities (Cij). Computer scientists call it an NP-complex problem [66], because it cannot be solved by an algorithm in polynomial time.Mathematically, TSP can be expressed as follows:Objective function:(1)min∑i,jcij⋅xij.Subject to:(2)∑ixij=1∀j≠i,∑jxij=1∀i≠j,u1=1,(3)2≤ui≤n∀i≠1,(4)ui−uj+1≤n−11−xij∀i≠1,∀j≠1,(5)ui≥0∀i,xij∈0,1∀i,j.Constraints (2)–(5) reach together to eliminate the subpaths in the solution.TSP can be used to solve different problems. Robotics can minimize the number of machine movements, where each movement has an associated cost [67]. Distribution logistics can find the minimum path from the point of origin, passing through known places and finally returning to the point of origin [68]. This happens in a dispatch system that minimizes fuel usage and ensures that the products are delivered as quickly as possible.


## 2.2. Honeybee Optimization Metaheuristics

Among the methods inspired by nature is the artificial bee colony (ABC) algorithm, which belongs to the swarm intelligence family. Within the bee colony methods, there are two distinct types of metaheuristics. The first approach is the bee collecting pollen algorithm (BCPA), where bees seek to obtain the maximum amount of nectar [69].MBO is based on the social behavior of bees in which queens, drones, larvae, and workers interact to solve overly complex combinatorial problems.


### 2.2.1. Description of the Metaheuristics

A bee colony can be formed in two different ways. The first method is called independent foundation, in which the colony is created by one or more queens that build a nest, produce eggs, and feed the larvae. Later, the queens only take charge of laying eggs, and the workers conduct the work in the hive. The second method is called swarming, in which one or more queens and workers create the colony from the original colony, and the division of labor begins immediately; therefore, the queens are only responsible for laying eggs. Colony formation by a single queen is called haplometrosis; otherwise, it is called pleometrosis. Additionally, if the hive has only one queen during its life cycle, it is monogynous, and if it has many queens, it is polygynous [70].Colonies are made up of queens, drones, larvae, and workers. Queens lay eggs, and male drones are haploid and can only affect the queen’s genotype via mutations. Fertilization occurs in the air, and the queen interbreeds 7 to 20 drones in one session. The drone’s genetic material is stored in the queen’s spermatheca, and a random mixture of the stored sperm fertilizes the queen’s eggs. The workers who care for the larvae also occasionally lay eggs.Our model randomly creates an initial population of queens, and its solutions are improved using workers. Each queen initiates a mating flight with initial energy and speed.The mating flight is a set of transitions where the queen moves in the air and mates (or does not mate) with a drone according to a probabilistic function. If mating occurs, the drone’s sperm is stored within the queen’s spermatheca. The flight ends when the queen’s energy has been exhausted or the spermatheca is at maximum capacity. After that, her energy and speed decreased.When the queen is in the nest, the drone’s genetic material is selected to create larvae cared for by the workers. The workers are represented by heuristics that seek to improve the larva (solution).Eventually, the best larvae become queens and replace the worst queens. The remaining queens and larvae are eliminated to begin a new flight cycle.The probability that queens and drones cross is represented by the annealing function, which depends on the queen’s fitness, the drone’s fitness, and the queen’s speed which is given as follows:(6)pq,d=mín1,e−lq,d/sq,where l (q, d) = dist(f(q), f(d)) corresponds to the absolute difference between the fitness of the queen and the drone, s(q) is the speed of the queen at that instant.The probability of mating is greater at the beginning of the flight since the queen’s energy is greater, and it is more likely the drone’s fitness is equivalent to the queen’s.As mentioned above, after each transition, speed and energy must be reduced. The speed must be adjusted according to the equation:(7)st+1=α⋅st,where α is a factor in the interval [0, 1].When energy is reduced, the following equation is used:(8)et+1=et−step,where step corresponds to the energy reduction, whose value is calculated according to the following equation:(9)step=0,5⋅e0M.M is the capacity of the spermatheca.Figure1 shows the five main stages of MBO: (1) mating flight in which the queen probabilistically selects the drones she will mate with and adds their sperm to her spermatheca; (2) creation of larvae by crossing between the genes of the queen and the drones; (3) using workers (heuristics) to improve the solution represented by the larvae; (4) update the worker’s fitness according to how well they improve the larvae; and (5) replace the worst queens with the best larvae.Figure 1
Genotype example.Abbass [31, 71] only considered a single queen and a single worker to solve a particular class of the satisfaction problem (3-SAT) and compared the results obtained with the WalkSAT and GSAT heuristics. MBO with the worker WalkSAT was more satisfying than WalkSAT alone, and so was MBO with the GSAT worker. The whole experiment was developed with parameters taken from biology.Teo and Abbass [72] changed the accepted trajectory between the drone and the queen. This time, the drone was accepted, and its sperm was added to the queen’s spermatheca if it was fitter than the queen or fulfilled the new annealing function which is given as follows:(10)pq,d=eΔf/sq,where Δ(f) corresponds to the difference between the fitness of the new trajectory and the queen and s(q) is the speed of the queen at that instant.Later, Abbass and Teo [36] changed the annealing function. Instead of using the difference between the queen’s fitness and the drone, they used the difference between the drone trajectories’ fitness to ensure that a trajectory was only accepted if it was better than the previous ones used. The new annealing function was determined by the following equation:(11)pq,d=e−ft−ft−1/sq,where (f(t) − f(t − 1)) corresponds to the difference in the fitness of the current trajectory and the previous trajectory and s(q) is the speed of the queen at that instant.The study addressed the 3-SAT problem using five workers: GSAT, random walk, probabilistic greedy, one-point crossover, and WalkSAT. The tests were conducted by taking each heuristic separately, using MBO with each worker, and using the workers together in MBO. The results showed that MBO behaves better when all the workers are used together and that the new modification [36] provides better results than the original version [31, 71] and the one proposed by [72].


### 2.2.2. State of the Art MBOs

Phu-Ang and Thammano [73] proposed a new MBO-based memetic algorithm to solve the flexible workshop scheduling problem. They introduced four new features to the standard MBO algorithm to move the search away from the local optimum. Their proposed algorithm maintained two populations: the female and male, and their performance was compared with the dispersion search algorithms using path linkage, hybrid differential evolution, hybrid harmony search algorithms, and hybrid GA. MBO outperformed others on large and complex problem cases.Wen et al. [74] proposed a modified HBMO (MHBMO) to integrate process planning and scheduling (IPPS) with uncertain processing time and due date based on the fuzzy set by designing an uncertainty measurement objective calculation method. Their method solves the IPPS of multiple objectives effectively.Palominos et al. [75] proposed a new MBO metaheuristic to deal with shop flow scheduling problems by designing different mating flight spaces and achieved excellent results over the 120 instances they tested.Solgi and Loáiciga [76] identified seven algorithms for solving continuous optimization problems and showed that ABC, bee evolution for genetic algorithms (BEGA), and MBO were the most efficient.Vakil-Baghmisheh and Salim [77] proposed a modified version of MFMHBO and compared it with ABC, QBE, and FMHBO in four reference functions for various variables up to 100. MFMHBO was faster than the others in most cases, especially for the Griewank and Schwefel functions, and increased precision and the number of variables. In general, their modified algorithm is comparable to the other algorithms.Prabhakar and Lee [78] studied different MBO techniques for transformation-based three-level character selection using wavelets for prostate cancer classification.Çelik and Ülker [79] tried to solve TSP using MBO. TSP had previously been solved using heuristic methods (GA and SA), but their MBO algorithm converged toward the optimal solution in fewer iterations.Nayak and Naik [80] presented a hybrid metaheuristic pi-sigma neural network (PSNN) based on honey bee mating that successfully solved the data mining classification problem. Their approach combined HBMO with PSNN and compared it with other techniques such as GA, differential evolution (DE), and PSO. Their approach was stable and reliable and provided better classification precision than others.Yin et al. [81] designed an accurate and stable MBO algorithm for programming a scorching machine to minimize the weighted sum of work completion times.Zanbouri and Jafari Navimipour [82] proposed creating cloud computing services where services were selected based on trust and quality (QoS). They used the media grouping algorithm k to reduce the search space and HBMO to achieve global optimization. Their method worked efficiently in terms of computational time, producing more reliable services, although it was not efficient for large-scale data sets.Arun and Vijay Kumar [83] proposed a modified MBO for a view selection algorithm (MBOVSA) to select views with the lowest total cost and information processing time in a large data warehouse. MBOVSA selected better quality views than HRUA, one of the important greed-based view selection algorithms.


## 2.2.1. Description of the Metaheuristics

A bee colony can be formed in two different ways. The first method is called independent foundation, in which the colony is created by one or more queens that build a nest, produce eggs, and feed the larvae. Later, the queens only take charge of laying eggs, and the workers conduct the work in the hive. The second method is called swarming, in which one or more queens and workers create the colony from the original colony, and the division of labor begins immediately; therefore, the queens are only responsible for laying eggs. Colony formation by a single queen is called haplometrosis; otherwise, it is called pleometrosis. Additionally, if the hive has only one queen during its life cycle, it is monogynous, and if it has many queens, it is polygynous [70].Colonies are made up of queens, drones, larvae, and workers. Queens lay eggs, and male drones are haploid and can only affect the queen’s genotype via mutations. Fertilization occurs in the air, and the queen interbreeds 7 to 20 drones in one session. The drone’s genetic material is stored in the queen’s spermatheca, and a random mixture of the stored sperm fertilizes the queen’s eggs. The workers who care for the larvae also occasionally lay eggs.Our model randomly creates an initial population of queens, and its solutions are improved using workers. Each queen initiates a mating flight with initial energy and speed.The mating flight is a set of transitions where the queen moves in the air and mates (or does not mate) with a drone according to a probabilistic function. If mating occurs, the drone’s sperm is stored within the queen’s spermatheca. The flight ends when the queen’s energy has been exhausted or the spermatheca is at maximum capacity. After that, her energy and speed decreased.When the queen is in the nest, the drone’s genetic material is selected to create larvae cared for by the workers. The workers are represented by heuristics that seek to improve the larva (solution).Eventually, the best larvae become queens and replace the worst queens. The remaining queens and larvae are eliminated to begin a new flight cycle.The probability that queens and drones cross is represented by the annealing function, which depends on the queen’s fitness, the drone’s fitness, and the queen’s speed which is given as follows:(6)pq,d=mín1,e−lq,d/sq,where l (q, d) = dist(f(q), f(d)) corresponds to the absolute difference between the fitness of the queen and the drone, s(q) is the speed of the queen at that instant.The probability of mating is greater at the beginning of the flight since the queen’s energy is greater, and it is more likely the drone’s fitness is equivalent to the queen’s.As mentioned above, after each transition, speed and energy must be reduced. The speed must be adjusted according to the equation:(7)st+1=α⋅st,where α is a factor in the interval [0, 1].When energy is reduced, the following equation is used:(8)et+1=et−step,where step corresponds to the energy reduction, whose value is calculated according to the following equation:(9)step=0,5⋅e0M.M is the capacity of the spermatheca.Figure1 shows the five main stages of MBO: (1) mating flight in which the queen probabilistically selects the drones she will mate with and adds their sperm to her spermatheca; (2) creation of larvae by crossing between the genes of the queen and the drones; (3) using workers (heuristics) to improve the solution represented by the larvae; (4) update the worker’s fitness according to how well they improve the larvae; and (5) replace the worst queens with the best larvae.Figure 1
Genotype example.Abbass [31, 71] only considered a single queen and a single worker to solve a particular class of the satisfaction problem (3-SAT) and compared the results obtained with the WalkSAT and GSAT heuristics. MBO with the worker WalkSAT was more satisfying than WalkSAT alone, and so was MBO with the GSAT worker. The whole experiment was developed with parameters taken from biology.Teo and Abbass [72] changed the accepted trajectory between the drone and the queen. This time, the drone was accepted, and its sperm was added to the queen’s spermatheca if it was fitter than the queen or fulfilled the new annealing function which is given as follows:(10)pq,d=eΔf/sq,where Δ(f) corresponds to the difference between the fitness of the new trajectory and the queen and s(q) is the speed of the queen at that instant.Later, Abbass and Teo [36] changed the annealing function. Instead of using the difference between the queen’s fitness and the drone, they used the difference between the drone trajectories’ fitness to ensure that a trajectory was only accepted if it was better than the previous ones used. The new annealing function was determined by the following equation:(11)pq,d=e−ft−ft−1/sq,where (f(t) − f(t − 1)) corresponds to the difference in the fitness of the current trajectory and the previous trajectory and s(q) is the speed of the queen at that instant.The study addressed the 3-SAT problem using five workers: GSAT, random walk, probabilistic greedy, one-point crossover, and WalkSAT. The tests were conducted by taking each heuristic separately, using MBO with each worker, and using the workers together in MBO. The results showed that MBO behaves better when all the workers are used together and that the new modification [36] provides better results than the original version [31, 71] and the one proposed by [72].


## 2.2.2. State of the Art MBOs

Phu-Ang and Thammano [73] proposed a new MBO-based memetic algorithm to solve the flexible workshop scheduling problem. They introduced four new features to the standard MBO algorithm to move the search away from the local optimum. Their proposed algorithm maintained two populations: the female and male, and their performance was compared with the dispersion search algorithms using path linkage, hybrid differential evolution, hybrid harmony search algorithms, and hybrid GA. MBO outperformed others on large and complex problem cases.Wen et al. [74] proposed a modified HBMO (MHBMO) to integrate process planning and scheduling (IPPS) with uncertain processing time and due date based on the fuzzy set by designing an uncertainty measurement objective calculation method. Their method solves the IPPS of multiple objectives effectively.Palominos et al. [75] proposed a new MBO metaheuristic to deal with shop flow scheduling problems by designing different mating flight spaces and achieved excellent results over the 120 instances they tested.Solgi and Loáiciga [76] identified seven algorithms for solving continuous optimization problems and showed that ABC, bee evolution for genetic algorithms (BEGA), and MBO were the most efficient.Vakil-Baghmisheh and Salim [77] proposed a modified version of MFMHBO and compared it with ABC, QBE, and FMHBO in four reference functions for various variables up to 100. MFMHBO was faster than the others in most cases, especially for the Griewank and Schwefel functions, and increased precision and the number of variables. In general, their modified algorithm is comparable to the other algorithms.Prabhakar and Lee [78] studied different MBO techniques for transformation-based three-level character selection using wavelets for prostate cancer classification.Çelik and Ülker [79] tried to solve TSP using MBO. TSP had previously been solved using heuristic methods (GA and SA), but their MBO algorithm converged toward the optimal solution in fewer iterations.Nayak and Naik [80] presented a hybrid metaheuristic pi-sigma neural network (PSNN) based on honey bee mating that successfully solved the data mining classification problem. Their approach combined HBMO with PSNN and compared it with other techniques such as GA, differential evolution (DE), and PSO. Their approach was stable and reliable and provided better classification precision than others.Yin et al. [81] designed an accurate and stable MBO algorithm for programming a scorching machine to minimize the weighted sum of work completion times.Zanbouri and Jafari Navimipour [82] proposed creating cloud computing services where services were selected based on trust and quality (QoS). They used the media grouping algorithm k to reduce the search space and HBMO to achieve global optimization. Their method worked efficiently in terms of computational time, producing more reliable services, although it was not efficient for large-scale data sets.Arun and Vijay Kumar [83] proposed a modified MBO for a view selection algorithm (MBOVSA) to select views with the lowest total cost and information processing time in a large data warehouse. MBOVSA selected better quality views than HRUA, one of the important greed-based view selection algorithms.


## 2.3. Experimental Settings

### 2.3.1. TSP Modelling using MBO

Our technique is based on [30]. We implemented Abbass’ proposal with some modifications in TSP. The MBO metaheuristic is described in Algorithm 1.Algorithm 1: MBO Metaheuristics.
Define the number of queensQ, Workers W, larvae B, mating flights G, and spermatheca sizeDefine the energy and speed of the queensAssign each worker with a different heuristicAssign the genotype of the queens randomlyChoose a worker randomly and apply its genotype to the queen’s genotypeWhile the number of mating flights<= GFor queen = 1 throughQAssign queen power and speedAssign step = 0.5 × energy/MAssign path randomlyAs long as energy> 0 and the spermatheca is not fullThe queen moves between different states and probabilistically chooses dronesIf a drone is selected, thenAdd your sperm to the queen’s spermathecaEndUpgrade energy and speed of the queenWhileEnd ForGenerate larva by crossing and mutationUse workers to improve larvaUpdate the fitness of the workersWhile the best of larvae is better than the worst of queensReplace the worst queen with the best larvaRemove the best larva from the larvae listEnd whileRemove all larvaEnd While


### 2.3.2. Defining Parameters

The proposed algorithm requires a series of parameters: (1) number of queen bees; (2) number of drones; (3) capacity of the queen’s spermatheca; (4) number of flights; (5) number of larvae; (6) initial energy; and (7) initial speed of the queen. Particular care was taken when choosing (7) (defined experimentally) because the speed influences the annealing function, which selects the genetic material of the drone. When (7) is low, it is more likely that a drone will not be chosen for a flight, and if it is high, it is more likely that all drones will be accepted. Additionally, since the fitness range can be broad in TSP, one should not select a fixed speed for all instances to avoid the above-mentioned problems. Therefore, (7) must be calculated considering the difference between the best and the worst queen when starting the algorithm. For parameter (6), values from previous studies published in [30] were used for the energy reduction and speed factors. Table 1 shows the parameter values used for the proposed MBO algorithm.Table 1
Algorithm parameters.
ParametersSource of obtainingValueNumber of queen beesExperimental3Drone quantityExperimental200Spermatheque capacityExperimental50Number of flightsExperimental100Number of larvaeExperimental30Initial energy[84]0.9Initial velocityExperimentalDifference between the best and the worst queen at the start of the algorithm.Energy reduction factor[30]step=0.5⋅e0/MSpeed reduction factor[30]α = 0.9


### 2.3.3. Representing the Solution

The genotype stereotype will be the same for queens, drones, and larvae and corresponds to a sequence of cities to be visited by the traveler. The genotype length is the number of cities that must be stopped by and varies according to the problem. Each gene contains the cities and the corresponding order. For example, if there is a three in position four, the city three will be after the fourth.Figure2 shows an example solution for a problem with six cities where the optimal order must be 1, 3, 6, 2, 5, and 4.Figure 2
HBMO algorithm [71].Each genotype has an associated fitness that corresponds to the tour distance.


### 2.3.4. Generating an Initial Population of Bees, Queens, and Drones

The initial generation of queens can be selected in three diverse ways. The first two queens are generated on the basis of two tour construction heuristics for TSP: the nearest neighbor [85] and the closest insertion [86]. The remaining queens and drones are randomly generated to introduce diversity [36].


### 2.3.5. Selection of Drone Genetic Material

Drone genetic material is selected by the MBO annealing function [30]. Drone (d) selection is more likely when its fitness is like the queen (q) and when the queen begins her flight.A random number between 0 and 1 is generated to select the drone. When this number is less than or equal to the one provided by the annealing function, the drone’s genetic material is added to the queen’s spermatheca (if the random number ≤pq,d, then add genetic material from the drone to the queen’s spermatheca).


### 2.3.6. Mating

After completing a flight, each queen returns to the nest to lay eggs. When all the queens have returned, the genetic materials used to generate larvae are selected. The number of larvae to be generated is an algorithm parameter, and each larva can come from any queen and any drone possessed by the chosen queen.Queens are chosen by fitness, i.e., fitter queens are more likely to be selected. Drone genetic material is picked up randomly, and larvae are generated by crossing.The crossing steps are as outlined below:(1)
Drones can only contribute half their genetic material because they are haploid [30, 87, 88]. Half of the drone’s genes tag like this to ensure.(2)
A random interval and segment (i, j), where i < j. i represents the beginning of the segment, and j represents the end.(3)
Each queen’s genes are checked, and they will be selected if they were not considered in the marked drone’s interval. The set of genes chosen by the queen generates a list.(4)
The first selectedi − 1 queen genes become part of the first i − 1 genes of the larva.(5)
The unlabeled genes in the marked drone segment (i, j) become part of the larva at the same position as the original gene.(6)
The queen’s gene that has not been considered goes back in the same position as the larva.(7)
From positionj + 1 the larva contains all other non-considered queen genes in the same order.Figure3 shows an example of the crossing process between a drone and a queen.Figure 3
Crossing.


### 2.3.7. Larva Mutation

The larva mutates after being formed by the following steps:(1)
Each gene is visited, and a random number between 0 and 1 is generated.(2)
The gene is marked if the random number is less than or equal to a mutation probability.(3)
Marked genes are mutated.We defined the mutation probability defined as 1% per [71].


### 2.3.8. Breeding Larvae Using Workers

We considered three types of workers corresponding to SA [89], path-based tabu search neural network (TSNN) [90], and swap local search (SLS) [91]. TSP created crossovers in our solutions, so we included an additional worker called the 2-opt [92] heuristic to avoid them. The assignment of this worker offers the solutions provided by the other three workers who were not chosen.


### 2.3.9. Updating the Queen Bee Population

The queen population is updated after a flight is finished. We use the same process as the original MBO algorithm (Algorithm2).Algorithm 2: Update the original queen population.
While the best of larvae is better than the worst of queensReplace the worst queen with the best larvaRemove the best larva from the larvae listEnd whileThe solutions with the fittest queens and larvae generate new queens. The remaining solutions are discarded but can be considered in the next iteration. A new flight begins after the queen population has been updated, unless the stop criterion is met.


### 2.3.10. Chaotic MBO Modeling

Chaos can be a useful tool for finding better solutions, but which MBO processes would benefit from chaos theory?The parts of the MBO metaheuristic (Algorithm1) must be considered, and a proposal to introduce chaos must be defined to answer this question. Chaos has successfully created initial solutions in other heuristics [93, 94].Drones are generated like queens, for which chaos better explores the solution space.In MBO, the workers improve the solutions. The effect of chaos cannot be evaluated if new chaotic heuristics replace MBO ones. Instead, chaos must be introduced into specific parts of the workers.We carried out different types of tests to measure the performance of chaos in the different MBO stages: (1) tests with traditional MBO; (2) tests with the use of chaos only in the generation of initial solutions (queens and drones); (3) tests with the use of chaos only in the workers; and (4) tests with the use of chaos in the generation of initial solutions and the workers.


### 2.3.11. Chaotic Generation of Initial Queen and Drone Populations

The initial populations of queen bees and the generation of the drone population are generated using chaos. One queen is born using the nearest neighbor heuristic, and another uses the closest insertion heuristic. There are two heuristics in both versions to ensure good solutions and improve the algorithm’s performance.The remaining queens are born using chaotic maps from a randomly initialized solution vector where each gene takes a value between 0 and 1 for each value within the logistic map solution as shown in the following equation:(12)Xn+1=rXn1−Xn,where r is the control parameter (in this case, r = 4), whose value ensures that the system will enter a state of chaos. r has values between 0 and 1, except for the fixed points of 0.25, 0.5, and 0.75.When a queen starts her flight, the first drone drops randomly and has a value between 0 and 1. The remaining drones are on the logistic map the same way as the queens.


### 2.3.12. Chaotic Breeding of Larvae using Workers

After generating the queens and larvae by crossing and mutations, the workers improve the solution quality.We use the same three workers with chaotic modifications (chaos simulated annealing (CSA), chaotic tabu search neural network based on a path (CTSNN), and chaotic swap local search (CSLS). Each worker performs a job that leads to better solutions.We included a chaotic version of the traditional MBO to avoid crossing TSPs in the solutions. This 2-opt heuristic uses the output of any chaotic worker mentioned above. The chaotic workers use the same fitness value at the beginning of the algorithm.


### 2.3.13. Chaos Simulated Annealing (CSA) Heuristic

The first chaotic worker follows the same idea of simulated annealing (SA) in the traditional MBO modeling but with chaotic alterations (CSA).CSA was created by Mingjun and Huan Wen [95] to take advantage of SA in combinatorial optimization problems and the improvements it brings to chaos theory.The problem:(13)Pminfx,xi∈ai,bi,i=1,2,…,n.CSA and SA have the same structure but differ in two ways. The first difference is when creating an initial solutionx0 = (x01, x02,…, x0n) using the logistic map to obtain zki, according to following equation:(14)x0i=ai+bi−ai⋅zki,i=1,…,n.The second difference is when selecting a neighboring solutionym = (xm1, xm2,…, xmn), where one of the current solution variables is random and it´s modified using the following logistic map again:(15)ymi=xmi+α⋅bi−ai⋅≤zkm,where α is a chaos weighting parameter.Algorithm3 shows the CSA heuristic proposed by [95].Algorithm 3: Pseudocode for CSA.
Chaotic initialization, giving an initial value toz0Generate the chaotic variableszki, i = 1, 2,…, n, using the logistic functionThe initial solutionx0 = (x01, x02,…,x0n) is produced by the formula x0i = ai + (bi − ai) × zki, i = 1,...,n.Initialize maximum temperatureTmax, minimum temperature Tmin and the number of iterations LmaxT =  Tmaxm = 0Set as the best solutionx∗=x0Set as the best fitnessf∗=fx∗While (T > Tmin)While (m ≤ Lmax)Generate a new solutionym = (xm1, xm2,…,xmn) based on the formula ymi = xmi + α × (bi − ai) × zkm, picking out i randomly and zkm is delivered by the logistics mapα = α × e-βΔE∗  = f(ym) − f∗ΔE = f(ym) − f(xm)If (ΔE∗ ≤ 0)x∗  =  ymf∗=fx∗If (ΔE ≤ 0)xm+1 =  ymf∗=fx∗If (ΔE > 0)Accept status with probabilitye−ΔE/Tm = m + 1End whileLmax  =  Lmax + d, m = 0T = δ × TEnd whileDeliver the best solution foundx∗ and best fitness f∗We worked on chaotic genotypes where each gene has a chaotic variable between 0 and 1. However, not all solutions have a chaotic genotype, such as the larvae born after applying the crossover and mutation operators, which generate a solution based on the traditional MBO. The genotype is a chaotic solution.Furthermore, chaotic MBO uses the same parameters as traditional MBO. Still, it is necessary to define the value of the chaos weighting (α), which we determined as one by experimentation.It is also essential to define the chaotic function used, which is determined by the original CSA study and corresponds to the logistic map:(16)zk+1=4zk1−zk.


### 2.3.14. Path-Based Chaotic Tabu Search Neural Network Heuristic

The second worker used for traditional MBO is path-based TSNN. Aono et al. [96] proposed a method based on the studies of [90], which incorporates chaotic noise into the equations.The new algorithm also relies on 2-opt moves to generate its TSP solutions. The new chaotic noise procedure is based on equations (17)–(19).(17)ξijt+1=βΔijt,(18)ζijt+1=−α∑d=0s−1krdxijt−d,(19)xijt+1=ξijt+1+ζaiajt+1+ζijt+1+γzijt+1,where Δij is the difference between the fitness of the current solution (D0) and the solution with a 2-opt movement between cities i and j (Dij): Δij (t) =  D0(t)  − Dij(t), β is the weight given to Δij(t), α is the weighting parameter of the taboo effect, kr is the decay parameter of the taboo effect, xij(t) is the output of neuron i, j in time t, and γ is the weight given to the chaotic noise zij.According to Aono et al. [96], the logistic map should generate zij in a chaotic form, as shown in the following equation:(20)zijt+1=3.828zijt1−zijt.The same parameters can be used for both the path-based TSNN and its chaotic version (CTSNN). The tests show different values ofγ in the interval [20, 30000], generating the best results for γ = 10000.


### 2.3.15. Chaotic Swap Local Search Heuristic

Like the chaotic heuristic described above, chaotic noise gets the worker to swap local search and creates the chaotic swap local search.To develop this new method, all swap movements (i, j), where (i, j) = (j, i) e i ≠ j, as in the following equation:(21)ξijt+1=Δijt+γzijt,where Δij is the difference between the values of the current solution’s objective function (D0) and the solution with a swap movement of cities i and j (Dij): Δij (t) =  D0(t)  − Dij(t) and γ is the weight given to the chaotic noise zij.The next movement selected corresponds to the one with the highest value ofξij, which is not necessarily the one that delivers the best fitness because of the chaotic variable, but the one with the swap movement between cities i and j. Finally, the solution with the best fitness is generated when all iterations are carried out.This worker uses the same parameters as the nonchaotic version, adding chaotic noise through the following equation:(22)zijt+1=3.828zijt1−zijt.We selected the chaos weighting experimentally, taking values ofγ in the interval [20, 100] and obtaining the best results for γ = 40.


### 2.3.16. Chaotic MBO Flow Chart

Figure4 shows a flow chart describing the chaotic MBO metaheuristic procedure described above to solve TPS.Figure 4
Chaotic MBO algorithm flow chart.


## 2.3.1. TSP Modelling using MBO

Our technique is based on [30]. We implemented Abbass’ proposal with some modifications in TSP. The MBO metaheuristic is described in Algorithm 1.Algorithm 1: MBO Metaheuristics.
Define the number of queensQ, Workers W, larvae B, mating flights G, and spermatheca sizeDefine the energy and speed of the queensAssign each worker with a different heuristicAssign the genotype of the queens randomlyChoose a worker randomly and apply its genotype to the queen’s genotypeWhile the number of mating flights<= GFor queen = 1 throughQAssign queen power and speedAssign step = 0.5 × energy/MAssign path randomlyAs long as energy> 0 and the spermatheca is not fullThe queen moves between different states and probabilistically chooses dronesIf a drone is selected, thenAdd your sperm to the queen’s spermathecaEndUpgrade energy and speed of the queenWhileEnd ForGenerate larva by crossing and mutationUse workers to improve larvaUpdate the fitness of the workersWhile the best of larvae is better than the worst of queensReplace the worst queen with the best larvaRemove the best larva from the larvae listEnd whileRemove all larvaEnd While


## 2.3.2. Defining Parameters

The proposed algorithm requires a series of parameters: (1) number of queen bees; (2) number of drones; (3) capacity of the queen’s spermatheca; (4) number of flights; (5) number of larvae; (6) initial energy; and (7) initial speed of the queen. Particular care was taken when choosing (7) (defined experimentally) because the speed influences the annealing function, which selects the genetic material of the drone. When (7) is low, it is more likely that a drone will not be chosen for a flight, and if it is high, it is more likely that all drones will be accepted. Additionally, since the fitness range can be broad in TSP, one should not select a fixed speed for all instances to avoid the above-mentioned problems. Therefore, (7) must be calculated considering the difference between the best and the worst queen when starting the algorithm. For parameter (6), values from previous studies published in [30] were used for the energy reduction and speed factors. Table 1 shows the parameter values used for the proposed MBO algorithm.Table 1
Algorithm parameters.
ParametersSource of obtainingValueNumber of queen beesExperimental3Drone quantityExperimental200Spermatheque capacityExperimental50Number of flightsExperimental100Number of larvaeExperimental30Initial energy[84]0.9Initial velocityExperimentalDifference between the best and the worst queen at the start of the algorithm.Energy reduction factor[30]step=0.5⋅e0/MSpeed reduction factor[30]α = 0.9


## 2.3.3. Representing the Solution

The genotype stereotype will be the same for queens, drones, and larvae and corresponds to a sequence of cities to be visited by the traveler. The genotype length is the number of cities that must be stopped by and varies according to the problem. Each gene contains the cities and the corresponding order. For example, if there is a three in position four, the city three will be after the fourth.Figure2 shows an example solution for a problem with six cities where the optimal order must be 1, 3, 6, 2, 5, and 4.Figure 2
HBMO algorithm [71].Each genotype has an associated fitness that corresponds to the tour distance.


## 2.3.4. Generating an Initial Population of Bees, Queens, and Drones

The initial generation of queens can be selected in three diverse ways. The first two queens are generated on the basis of two tour construction heuristics for TSP: the nearest neighbor [85] and the closest insertion [86]. The remaining queens and drones are randomly generated to introduce diversity [36].


## 2.3.5. Selection of Drone Genetic Material

Drone genetic material is selected by the MBO annealing function [30]. Drone (d) selection is more likely when its fitness is like the queen (q) and when the queen begins her flight.A random number between 0 and 1 is generated to select the drone. When this number is less than or equal to the one provided by the annealing function, the drone’s genetic material is added to the queen’s spermatheca (if the random number ≤pq,d, then add genetic material from the drone to the queen’s spermatheca).


## 2.3.6. Mating

After completing a flight, each queen returns to the nest to lay eggs. When all the queens have returned, the genetic materials used to generate larvae are selected. The number of larvae to be generated is an algorithm parameter, and each larva can come from any queen and any drone possessed by the chosen queen.Queens are chosen by fitness, i.e., fitter queens are more likely to be selected. Drone genetic material is picked up randomly, and larvae are generated by crossing.The crossing steps are as outlined below:(1)
Drones can only contribute half their genetic material because they are haploid [30, 87, 88]. Half of the drone’s genes tag like this to ensure.(2)
A random interval and segment (i, j), where i < j. i represents the beginning of the segment, and j represents the end.(3)
Each queen’s genes are checked, and they will be selected if they were not considered in the marked drone’s interval. The set of genes chosen by the queen generates a list.(4)
The first selectedi − 1 queen genes become part of the first i − 1 genes of the larva.(5)
The unlabeled genes in the marked drone segment (i, j) become part of the larva at the same position as the original gene.(6)
The queen’s gene that has not been considered goes back in the same position as the larva.(7)
From positionj + 1 the larva contains all other non-considered queen genes in the same order.Figure3 shows an example of the crossing process between a drone and a queen.Figure 3
Crossing.


## 2.3.7. Larva Mutation

The larva mutates after being formed by the following steps:(1)
Each gene is visited, and a random number between 0 and 1 is generated.(2)
The gene is marked if the random number is less than or equal to a mutation probability.(3)
Marked genes are mutated.We defined the mutation probability defined as 1% per [71].


## 2.3.8. Breeding Larvae Using Workers

We considered three types of workers corresponding to SA [89], path-based tabu search neural network (TSNN) [90], and swap local search (SLS) [91]. TSP created crossovers in our solutions, so we included an additional worker called the 2-opt [92] heuristic to avoid them. The assignment of this worker offers the solutions provided by the other three workers who were not chosen.


## 2.3.9. Updating the Queen Bee Population

The queen population is updated after a flight is finished. We use the same process as the original MBO algorithm (Algorithm2).Algorithm 2: Update the original queen population.
While the best of larvae is better than the worst of queensReplace the worst queen with the best larvaRemove the best larva from the larvae listEnd whileThe solutions with the fittest queens and larvae generate new queens. The remaining solutions are discarded but can be considered in the next iteration. A new flight begins after the queen population has been updated, unless the stop criterion is met.


## 2.3.10. Chaotic MBO Modeling

Chaos can be a useful tool for finding better solutions, but which MBO processes would benefit from chaos theory?The parts of the MBO metaheuristic (Algorithm1) must be considered, and a proposal to introduce chaos must be defined to answer this question. Chaos has successfully created initial solutions in other heuristics [93, 94].Drones are generated like queens, for which chaos better explores the solution space.In MBO, the workers improve the solutions. The effect of chaos cannot be evaluated if new chaotic heuristics replace MBO ones. Instead, chaos must be introduced into specific parts of the workers.We carried out different types of tests to measure the performance of chaos in the different MBO stages: (1) tests with traditional MBO; (2) tests with the use of chaos only in the generation of initial solutions (queens and drones); (3) tests with the use of chaos only in the workers; and (4) tests with the use of chaos in the generation of initial solutions and the workers.


## 2.3.11. Chaotic Generation of Initial Queen and Drone Populations

The initial populations of queen bees and the generation of the drone population are generated using chaos. One queen is born using the nearest neighbor heuristic, and another uses the closest insertion heuristic. There are two heuristics in both versions to ensure good solutions and improve the algorithm’s performance.The remaining queens are born using chaotic maps from a randomly initialized solution vector where each gene takes a value between 0 and 1 for each value within the logistic map solution as shown in the following equation:(12)Xn+1=rXn1−Xn,where r is the control parameter (in this case, r = 4), whose value ensures that the system will enter a state of chaos. r has values between 0 and 1, except for the fixed points of 0.25, 0.5, and 0.75.When a queen starts her flight, the first drone drops randomly and has a value between 0 and 1. The remaining drones are on the logistic map the same way as the queens.


## 2.3.12. Chaotic Breeding of Larvae using Workers

After generating the queens and larvae by crossing and mutations, the workers improve the solution quality.We use the same three workers with chaotic modifications (chaos simulated annealing (CSA), chaotic tabu search neural network based on a path (CTSNN), and chaotic swap local search (CSLS). Each worker performs a job that leads to better solutions.We included a chaotic version of the traditional MBO to avoid crossing TSPs in the solutions. This 2-opt heuristic uses the output of any chaotic worker mentioned above. The chaotic workers use the same fitness value at the beginning of the algorithm.


## 2.3.13. Chaos Simulated Annealing (CSA) Heuristic

The first chaotic worker follows the same idea of simulated annealing (SA) in the traditional MBO modeling but with chaotic alterations (CSA).CSA was created by Mingjun and Huan Wen [95] to take advantage of SA in combinatorial optimization problems and the improvements it brings to chaos theory.The problem:(13)Pminfx,xi∈ai,bi,i=1,2,…,n.CSA and SA have the same structure but differ in two ways. The first difference is when creating an initial solutionx0 = (x01, x02,…, x0n) using the logistic map to obtain zki, according to following equation:(14)x0i=ai+bi−ai⋅zki,i=1,…,n.The second difference is when selecting a neighboring solutionym = (xm1, xm2,…, xmn), where one of the current solution variables is random and it´s modified using the following logistic map again:(15)ymi=xmi+α⋅bi−ai⋅≤zkm,where α is a chaos weighting parameter.Algorithm3 shows the CSA heuristic proposed by [95].Algorithm 3: Pseudocode for CSA.
Chaotic initialization, giving an initial value toz0Generate the chaotic variableszki, i = 1, 2,…, n, using the logistic functionThe initial solutionx0 = (x01, x02,…,x0n) is produced by the formula x0i = ai + (bi − ai) × zki, i = 1,...,n.Initialize maximum temperatureTmax, minimum temperature Tmin and the number of iterations LmaxT =  Tmaxm = 0Set as the best solutionx∗=x0Set as the best fitnessf∗=fx∗While (T > Tmin)While (m ≤ Lmax)Generate a new solutionym = (xm1, xm2,…,xmn) based on the formula ymi = xmi + α × (bi − ai) × zkm, picking out i randomly and zkm is delivered by the logistics mapα = α × e-βΔE∗  = f(ym) − f∗ΔE = f(ym) − f(xm)If (ΔE∗ ≤ 0)x∗  =  ymf∗=fx∗If (ΔE ≤ 0)xm+1 =  ymf∗=fx∗If (ΔE > 0)Accept status with probabilitye−ΔE/Tm = m + 1End whileLmax  =  Lmax + d, m = 0T = δ × TEnd whileDeliver the best solution foundx∗ and best fitness f∗We worked on chaotic genotypes where each gene has a chaotic variable between 0 and 1. However, not all solutions have a chaotic genotype, such as the larvae born after applying the crossover and mutation operators, which generate a solution based on the traditional MBO. The genotype is a chaotic solution.Furthermore, chaotic MBO uses the same parameters as traditional MBO. Still, it is necessary to define the value of the chaos weighting (α), which we determined as one by experimentation.It is also essential to define the chaotic function used, which is determined by the original CSA study and corresponds to the logistic map:(16)zk+1=4zk1−zk.


## 2.3.14. Path-Based Chaotic Tabu Search Neural Network Heuristic

The second worker used for traditional MBO is path-based TSNN. Aono et al. [96] proposed a method based on the studies of [90], which incorporates chaotic noise into the equations.The new algorithm also relies on 2-opt moves to generate its TSP solutions. The new chaotic noise procedure is based on equations (17)–(19).(17)ξijt+1=βΔijt,(18)ζijt+1=−α∑d=0s−1krdxijt−d,(19)xijt+1=ξijt+1+ζaiajt+1+ζijt+1+γzijt+1,where Δij is the difference between the fitness of the current solution (D0) and the solution with a 2-opt movement between cities i and j (Dij): Δij (t) =  D0(t)  − Dij(t), β is the weight given to Δij(t), α is the weighting parameter of the taboo effect, kr is the decay parameter of the taboo effect, xij(t) is the output of neuron i, j in time t, and γ is the weight given to the chaotic noise zij.According to Aono et al. [96], the logistic map should generate zij in a chaotic form, as shown in the following equation:(20)zijt+1=3.828zijt1−zijt.The same parameters can be used for both the path-based TSNN and its chaotic version (CTSNN). The tests show different values ofγ in the interval [20, 30000], generating the best results for γ = 10000.


## 2.3.15. Chaotic Swap Local Search Heuristic

Like the chaotic heuristic described above, chaotic noise gets the worker to swap local search and creates the chaotic swap local search.To develop this new method, all swap movements (i, j), where (i, j) = (j, i) e i ≠ j, as in the following equation:(21)ξijt+1=Δijt+γzijt,where Δij is the difference between the values of the current solution’s objective function (D0) and the solution with a swap movement of cities i and j (Dij): Δij (t) =  D0(t)  − Dij(t) and γ is the weight given to the chaotic noise zij.The next movement selected corresponds to the one with the highest value ofξij, which is not necessarily the one that delivers the best fitness because of the chaotic variable, but the one with the swap movement between cities i and j. Finally, the solution with the best fitness is generated when all iterations are carried out.This worker uses the same parameters as the nonchaotic version, adding chaotic noise through the following equation:(22)zijt+1=3.828zijt1−zijt.We selected the chaos weighting experimentally, taking values ofγ in the interval [20, 100] and obtaining the best results for γ = 40.


## 2.3.16. Chaotic MBO Flow Chart

Figure4 shows a flow chart describing the chaotic MBO metaheuristic procedure described above to solve TPS.Figure 4
Chaotic MBO algorithm flow chart.


## 3. Results

We carried out the computational experiment in the programming language C and extracted the used instances from the TSPLIB database [97] (Table 2). We solved each instance using the traditional MBO (MBO1) using a common pseudorandom number generator (PRNG) in the evolutionary operators and three chaotic MBO variants: (1) MBO with chaos in the initial solution generation (MBO2); (2) MBO with chaos in the workers (MBO3); and (3) MBO with chaos in the initial solution generation and the workers (MBO4). We executed each instance five times under the parameters already described in the previous points to find an average error for each one.Table 2
Description of the instances.
InstanceNumber of citiesKroA100100KroB100100KroC100100KroD100100KroE100100KroA150150KroB150150d198198KroA200200KroB200200a280280pr299299lin318318rd400400fl417417pcb442442d493492d657657rat783783fl15771577The solution quality is in equation (7), which is the percentage error between the average value of the solutions found (PSE) for the pure MBO and the three chaotic variants with the best-known solutions (BKS) of the TSPLIB as shown in the following equation:(23)Percent Error=PSE−BKSBKS∗100.The results for the four metaheuristics are shown in Table3. MBO3 and MBO4 provided the best results. However, ANOVA and Tukey statistical analyzes (with a 95% confidence interval) revealed that 30% of the instances were significantly different across the four methods applied.Table 3
Metaheuristic results MBO1, MBO2, MBO3, MBO4.
InstanceValue optimumMBO 1Error %MBO 2Error %MBO 3Error %MBO 4Error %KroA1002128221285.40.0221285.40.0221285.40.0221285.40.02KroB1002214122150.70.0422139.1−0.0122139.1−0.0122139.1−0.01KroC1002074920750.80.0120750.80.0120750.80.0120750.80.01KroD1002129421297.30.0221294.30.0021294.30.0021294.30.00KroE1002206822113.70.2122068.80.0022068.80.0022086.50.08KroA1502652426606.90.3126601.30.2926535.10.0426535.20.04KroB1502613026155.80.1026148.00.0726137.80.0326133.60.01d1981578015848.80.4415846.60.4215823.70.2815821.80.27KroA2002936829452.20.2929515.00.5029433.90.2229432.50.22KroB2002943729583.30.5029545.80.3729539.60.3529516.40.27a28025792603.00.932603.50.952588.50.372588.50.37pr2994819148578.20.8048467.80.5748409.50.4548389.60.41lin3184202942649.61.4842522.41.1742513.11.1542331.10.72rd4001528115637.92.3415651.02.4215582.81.9815590.52.03fl4171186111961.70.8511962.20.8511931.10.5911940.00.67pcb4425077951748.41.9151953.02.3151586.91.5951657.01.73d4933500235689.11.9635670.01.9135478.21.3635589.91.68d6574891250430.23.1050301.22.8450298.82.8450334.12.91rat78388069194.04.419186.84.329172.14.169236.54.89fl15772224922635.81.7422696.02.0122732.22.1722653.41.82Average error1.071.050.880.91Additionally, we found the chaotic modifications improved convergence speed (Figure5). MBO3 and MBO4 converged faster than MBO1.Figure 5
Convergence curve for instance d198.


### 3.1. Discussion of the Results

The calculated percentage error of pure MBO was 1.07% when it was compared with the results of 20 instances with the best values found in the literature [98]. When generating solutions, MBO with chaos in the workers and MBO with chaos in both have percentage errors of 1.05%, 0.88%, and 0.91%, respectively.In the pure MBO, the TSNN worker had the highest efficiency based on routes, followed by SA and later SLS. The CTSNN worker-generated the best results in the versions with chaos based on routes, then CSLS and CSA.At first sight, chaos improves solution quality when analyzed. The ANOVA analysis revealed six instances of significant differences between the four methods. The six TSP instances compared with the TSPLIB library are kroE100, d198, a280, pr299, fl417, and pcb442. Twenty experiments per metaheuristic were performed for each one.Following the results of the ANOVA test, the Tukey nonparametric evaluation was used to compare the two groups of metaheuristics. The results of the Tukey nonparametric test in five out of six instances show a statistical difference between metaheuristics with chaos and metaheuristics with MBO1. The results of the ANOVA test and the Tukey nonparametric test are shown in Tables4–15.Table 4
ANOVA for the kroE100 instance.
Origin of variationsSum of squaresDegrees of freedomAverage of the squaresFP valueCritical value forFBetween groups6772,796332257,598776,825129970,003572223,23887152Within groups5292,4384516330,777403———Total12065,234719————Table 5
Tukey’s test for the kroE100 instance.
TypeNSubset for alpha = .0512MBO2522068,8—MBO3522068,8—MBO4522086,5222086,52MBO15—22113,74p-value0,437085730,12380896Table 6
ANOVA for the d198 instance.
Origin of variationsSum of squaresDegrees of freedomAverage of the squaresFP valueCritical value forFBetween groups3136,3790131045,459677,455208140,002421983,23887152Within groups2243,7139816140,232124———Total5380,0929919————Table 7
Tukey’s test for the d198 instance.
TypeNSubset for alpha = .0512MBO4515821,84—MBO3515823,74—MBO25—15846,64MBO15—15848,82P value0,994036980,99107353Table 8
ANOVA for the a280 instance.
Origin of variationsSum of squaresDegrees of freedomAverage of the squaresFP valueCritical value forFBetween groups1081,40573360,46856739,09807871,3599E − 073,23887152Within groups147,513567169,21959795———Total1228,9192719————Table 9
Tukey’s test for the a280 instance.
TypeNSubset for alpha = 0.0512MBO352588,52—MBO452588,54—MBO15—2602,98MBO25—2603,46pvalue0,999999570,99426841Table 10
ANOVA for the pr299 instance.
Origin of variationsSum of squaresDegrees of freedomAverage of the squaresFP valueCritical value forFBetween groups107695,484335898,49473,59386040,037031613,23887152Within groups159821,432169988,83951———Total267516,91619————Table 11
Tukey’s test for the pr299 instance.
TypeNSubset for alpha = 0.0512MBO4548389,56—MBO3548409,548409,5MBO2548467,7848467,78MBO15—48578,2p value0,613333540,07194982Table 12
ANOVA for the fl417 instance.
Origin of variationsSum of squaresDegrees of freedomAverage of the squaresFP valueCritical value forFBetween groups3678,4089131226,136311,76115910,000254723,23887152Within groups1668,0482616104,253016———Total5346,4571719————Table 13
Tukey’s test for the fl417 instance.
TypeNSubset for alpha = 0.0512MBO3511931,1—MBO4511940,04—MBO15—11961,74MBO25—11962,22p value0,527135160,99984659Table 14
ANOVA for the pcb442 instance.
Origin of variationsSum of squaresDegrees of freedomAverage of the squaresFP valueCritical value forFBetween groups378562,5223126187,5073,436943770,042235583,23887152Within groups587440,5441636715,034———Total966003,06619————Table 15
Tukey’s test for the pcb442 instance.
TypeNSubset for alpha = 0.0512MBO3551586,96—MBO455165751657MBO1551748,4451748,44MBO25—51953,04p value0,556744310,10864054The results show that chaotic functions applied to the MBO metaheuristic significantly improve solutions regarding the TSP. All the solutions are equal to or greater than those obtained with pure MBO.


## 3.1. Discussion of the Results

The calculated percentage error of pure MBO was 1.07% when it was compared with the results of 20 instances with the best values found in the literature [98]. When generating solutions, MBO with chaos in the workers and MBO with chaos in both have percentage errors of 1.05%, 0.88%, and 0.91%, respectively.In the pure MBO, the TSNN worker had the highest efficiency based on routes, followed by SA and later SLS. The CTSNN worker-generated the best results in the versions with chaos based on routes, then CSLS and CSA.At first sight, chaos improves solution quality when analyzed. The ANOVA analysis revealed six instances of significant differences between the four methods. The six TSP instances compared with the TSPLIB library are kroE100, d198, a280, pr299, fl417, and pcb442. Twenty experiments per metaheuristic were performed for each one.Following the results of the ANOVA test, the Tukey nonparametric evaluation was used to compare the two groups of metaheuristics. The results of the Tukey nonparametric test in five out of six instances show a statistical difference between metaheuristics with chaos and metaheuristics with MBO1. The results of the ANOVA test and the Tukey nonparametric test are shown in Tables4–15.Table 4
ANOVA for the kroE100 instance.
Origin of variationsSum of squaresDegrees of freedomAverage of the squaresFP valueCritical value forFBetween groups6772,796332257,598776,825129970,003572223,23887152Within groups5292,4384516330,777403———Total12065,234719————Table 5
Tukey’s test for the kroE100 instance.
TypeNSubset for alpha = .0512MBO2522068,8—MBO3522068,8—MBO4522086,5222086,52MBO15—22113,74p-value0,437085730,12380896Table 6
ANOVA for the d198 instance.
Origin of variationsSum of squaresDegrees of freedomAverage of the squaresFP valueCritical value forFBetween groups3136,3790131045,459677,455208140,002421983,23887152Within groups2243,7139816140,232124———Total5380,0929919————Table 7
Tukey’s test for the d198 instance.
TypeNSubset for alpha = .0512MBO4515821,84—MBO3515823,74—MBO25—15846,64MBO15—15848,82P value0,994036980,99107353Table 8
ANOVA for the a280 instance.
Origin of variationsSum of squaresDegrees of freedomAverage of the squaresFP valueCritical value forFBetween groups1081,40573360,46856739,09807871,3599E − 073,23887152Within groups147,513567169,21959795———Total1228,9192719————Table 9
Tukey’s test for the a280 instance.
TypeNSubset for alpha = 0.0512MBO352588,52—MBO452588,54—MBO15—2602,98MBO25—2603,46pvalue0,999999570,99426841Table 10
ANOVA for the pr299 instance.
Origin of variationsSum of squaresDegrees of freedomAverage of the squaresFP valueCritical value forFBetween groups107695,484335898,49473,59386040,037031613,23887152Within groups159821,432169988,83951———Total267516,91619————Table 11
Tukey’s test for the pr299 instance.
TypeNSubset for alpha = 0.0512MBO4548389,56—MBO3548409,548409,5MBO2548467,7848467,78MBO15—48578,2p value0,613333540,07194982Table 12
ANOVA for the fl417 instance.
Origin of variationsSum of squaresDegrees of freedomAverage of the squaresFP valueCritical value forFBetween groups3678,4089131226,136311,76115910,000254723,23887152Within groups1668,0482616104,253016———Total5346,4571719————Table 13
Tukey’s test for the fl417 instance.
TypeNSubset for alpha = 0.0512MBO3511931,1—MBO4511940,04—MBO15—11961,74MBO25—11962,22p value0,527135160,99984659Table 14
ANOVA for the pcb442 instance.
Origin of variationsSum of squaresDegrees of freedomAverage of the squaresFP valueCritical value forFBetween groups378562,5223126187,5073,436943770,042235583,23887152Within groups587440,5441636715,034———Total966003,06619————Table 15
Tukey’s test for the pcb442 instance.
TypeNSubset for alpha = 0.0512MBO3551586,96—MBO455165751657MBO1551748,4451748,44MBO25—51953,04p value0,556744310,10864054The results show that chaotic functions applied to the MBO metaheuristic significantly improve solutions regarding the TSP. All the solutions are equal to or greater than those obtained with pure MBO.


## 4. Conclusions

The research and analysis of new methods for solving difficult optimization problems are of great significance for this type of concern. When introduced, chaos theory in MBO metaheuristics, specifically when generating solutions and in the workers, created four algorithms. The first considers traditional MBO without chaos; the second only uses chaos when generating solutions; the third uses chaos in the workers; and the fourth uses chaos when generating solutions and in the workers. In this study, TSP was used to measure the performance of each method.This research demonstrates the ability to improve performance using chaos theory. However, it is impossible to extrapolate the results obtained with this chaotic MBO heuristic to all engineering problems. The performance of each heuristic is closely related to the solution area of the problem; for the TSP problem, the results show better performance by introducing chaos to the general MBO model. The technique improved the global search capability of the MBO3 heuristic with a chaotic modification of the worker bees by concentrating on the highest number of high-performance solutions.Future work should focus on connecting MBO with chaos theory and using it in other stages of MBO, such as the mutation that increases larvae, the drone quality of life before mating, and new chaotic heuristics in the workers. New versions can be used by the technique using parallel distributions or by hybridizing with different chaotic mappings to improve performance.



---

*Source: 8903005-2022-10-11.xml*