# On the Convergence of Biogeography-Based Optimization for Binary Problems

Authors: Haiping Ma; Dan Simon; Minrui Fei
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

Biogeography-based optimization (BBO) is an evolutionary algorithm inspired by biogeography, which is the study of the migration of species between habitats. A finite Markov chain model of BBO for binary problems was derived in earlier work, and some significant theoretical results were obtained. This paper analyzes the convergence properties of BBO on binary problems based on the previously derived BBO Markov chain model. Analysis reveals that BBO with only migration and mutation never converges to the global optimum. However, BBO with elitism, which maintains the best candidate in the population from one generation to the next, converges to the global optimum. In spite of previously published differences between genetic algorithms (GAs) and BBO, this paper shows that the convergence properties of BBO are similar to those of the canonical GA. In addition, the convergence rate estimate of BBO with elitism is obtained in this paper and is confirmed by simulations for some simple representative problems.


---

## Body

## 1. Introduction

Mathematical models of biogeography describe the immigration and emigration of species between habitats. Biogeography-based optimization (BBO) was first presented in 2008 [1] and is an extension of biogeography theory to evolutionary algorithms (EAs). BBO is modeled after the immigration and emigration of species between habitats. One distinctive feature of BBO is that it uses the fitness of each candidate solution to determine its immigration and emigration rate. The immigration rate determines how likely a candidate solution is to change its decision variables, and the emigration rate determines how likely a candidate solution is to share its decision variables with other candidate solutions. Specifically, a candidate solution’s emigration rate increases with fitness, and its immigration rate decreases with fitness.Although BBO is a relatively new EA, it has demonstrated good performance on various unconstrained and constrained benchmark functions [2–5] and on real-world optimization problems such as sensor selection [1], economic load dispatch [6], robot controller tuning [7], satellite image classification [8], and power system optimization [9]. In addition, Markov models have been derived for BBO on binary problems [10, 11]. Reference [12] discusses the conceptual, algorithmic, and performance differences between BBO and GAs using both Markov model comparisons and benchmark simulation results. These simulation and theoretical results confirm that BBO is a competitive evolutionary algorithm. But until now there have not been any theoretical results concerning its convergence properties.We say that an optimization algorithm converges to the global optimum if the value of at least one of its candidate solutions, in the limit as the generation count approaches infinity, is equal to the global optimum of the optimization problem. Several mathematical tools have been used to analyze EA convergence [13–17]. Recent work includes the analysis of EA convergence using Markov’s inequality, Chebyshev bounds, Chernoff bounds, and martingales for minimum spanning tree problems, maximum matching problems, scheduling problems, shortest path problems, Eulerian cycle problems, multiobjective problems, and others [18–20].Markov chain models are still some of the most frequently used methods for the analysis of EA convergence. They have been widely used in a variety of EAs, including genetic algorithms (GAs) [21–25] and simulated annealing [4, 26], to prove probabilistic convergence to the global optimum. A Markov chain is a random process which has a discrete set of possible states si(i=1,2,…,T). The probability that the system transitions from state si to si is given by pij, which is called a transition probability. The T×T matrix P = [pij] is called the transition matrix, where T is the total number of possible population distributions. A population distribution is a specific multiset of individuals with a cardinality that is equal to the population size. A Markov state in [11] represents a BBO population distribution. Each state represents a particular population distribution, that is, how many individuals at each point of the search space there are in the population. Probability pij is the probability that the population transitions from the ith population distribution to the jth population distribution in one generation.This paper analyzes the global convergence properties of BBO as applied to optimization problems with binary search spaces, based on a previously derived BBO Markov model, and obtains the convergence rate estimate using homogeneous finite Markov chain properties. Section2 gives a brief review of BBO and its Markov transition probabilities. Section 3 gives some basic definitions, obtains BBOs convergence properties, and obtains the convergence rate estimate. Section 4 confirms the theory using simple numerical simulation. The convergence properties and convergence rates derived here are not surprising in view of previous EA convergence results, but this paper represents the first time that such results have been formalized for BBO. Finally, Section 5 presents some concluding remarks and directions for future work.


## 2. Biogeography-Based Optimization (BBO)

This section presents a review of the biogeography-based optimization (BBO) algorithm with migration and mutation (Section2.1) and then provides a review of the Markov transition probability of BBO populations (Section 2.2).


### 2.1. Overview of Biogeography-Based Optimization

This section provides an overview of BBO. The review in this section is very general because it applies to optimization problems with either real domains, integer domains, binary domains, or combinations thereof.BBO is a new optimization approach, inspired by biogeography theory, to solve general optimization problems. A biogeography habitat corresponds to acandidate solution to the optimization problem. A multiset of biogeography habitats corresponds to apopulation of candidate solutions. Habitat suitability index (HSI) in biogeography corresponds to the goodness of a candidate solution, which is also calledfitnessin standard EA notation. Like other EAs [27], BBO probabilistically shares information between candidate solutions to improve candidate solution fitness. In BBO, each candidate solution is comprised of a set of features, which are also called independent variables or decision variables in the optimization literature. Note that the decision variables can be taken from sets of real numbers, integers, binary numbers, or combinations thereof, depending on the problem. Each candidate solution immigrates features from other candidate solutions based on its immigration rate and emigrates features to other candidate solutions based on its emigration rate. BBO consists of two main steps: migration and mutation.Migration. Migration is a probabilistic operator that is intended to improve a candidate solution yk. For each feature of a given candidate solution yk, the candidate solution’s immigration rate λk is used to probabilistically decide whether or not to immigrate. If immigration occurs, then the emigrating candidate solution yj is probabilistically chosen based on the emigration rate μj. Migration is written as follows:
(1)yk(r)⟵yj(r),
where r is a candidate solution feature. In BBO, each candidate solution yk has its own immigration rate λk and emigration rate μk. A good candidate solution has relatively high μ and low λ, while the converse is true for a poor candidate solution. According to [1], the immigration rate λk and emigration rate μk of the candidate solution yk are based on linear functions and are calculated as
(2)λk=1-fitness(yk)μk=fitness(yk),
where fitness denotes candidate solution fitness value, which is normalized to the range [0, 1]. The probabilities of immigrating to yk and of emigrating from yk are calculated as
(3)Pr(immigrationtoyk)=λkPr(emigrationfromyk)=μk∑j=1Nμj,
where N is the population size.Mutation. Mutation is a probabilistic operator that randomly modifies a candidate solution feature. The purpose of mutation is to increase diversity among the population, just as in other EAs. Mutation of the ith candidate solution is implemented as shown in Algorithm 1.Algorithm 1
F
or each candidate solution feature
r  inyi
Ifrand(0,1)  <pm
yi(r)←  rand(Lr,Ur)
End if
Next
rIn the mutation logic (Algorithm1), rand (a,b) is a uniformly distributed random number between a and b, pm is the mutation rate, and Lr and Ur are the lower and upper search bounds of the rth independent variable. The above logic mutates each independent variable with a probability of pm. If mutation occurs for a given independent variable, then that independent variable is replaced with a random number within its search domain.A description of one generation of BBO is given in Algorithm2.Algorithm 2:One generation of the BBO algorithm. y is the entire population of candidate solutions, yk is the kth candidate solution, and yk(r) is the rth feature of yk.
For eachyk, define emigration rate μk proportional to fitness of yk, μk∈[0,1]
For eachyk, define immigration rate λk=1-μk
For eachyi
For each candidate solution featurer
Use λi to probabilistically decide whether to immigrate to yi (see (2) and (3))
If immigrating then
Use{μk} to probabilistically select the emigrating solution yk (see (2) and (3))
yi(r)←yk(r)
End if
Next candidate solution feature
Probabilistically decide whether to mutateyi
Next solution


### 2.2. Transition Probability of Biogeography-Based Optimization

Two main steps in BBO are significant: migration and mutation. The transition probability of one generation includes migration probability and mutation probability. Consider a problem whose search space{x1,…,xn} is binary. The search space is the set of all bit strings xi consisting of q bits each. Therefore, the cardinality of the search space is n=2q. Suppose that BBO is currently in the tth generation. Based on the previously derived Markov chain model for BBO [11], the probability that migration results in the kth candidate solution yk are equal to xi at generation t+1 is given by
(4)Pr(yk,t+1=xi)=∏r=1q[(1-λk)10(yk(r)-xi(r))+λk∑j∈ςi(r)vjμj∑j=1nvjμj],
where10 is the indicator function on the set 0 (i.e.,10(a)=1 if a = 0, and10(a) = 0 if a≠0), r denotes the index of the candidate solution feature (i.e., the bit number), λk denotes the immigration rate of candidate solution yk, μj denotes the emigration rate of candidate solution xj, and vj denotes the number of xj individuals in the population. The notation ςi(r) in (4) denotes the set of search space indices j such that the rth bit of xj is equal to the rth bit of xi. That is, ςi(r) = {j:xj(r)=xi(r)}. Note that the first term in the product on the right side of (4) denotes the probability that yk,t+1(r)=xi(r) if immigration of the rth candidate solution feature does not occur, and the second term denotes the probability if immigration of the rth candidate solution feature does occur.Example 1.
To clarify the notations in (4), an example is presented. We use the notation vj to denote the number of xj individuals in the population. Suppose we have a two-bit problem (q=2, n=4) with a population size N=3. The search space consists of the bit strings x={x1, x2, x3, x4} = {00, 01, 10, 11}. Suppose that the individuals in the current population are y={x2, x4, x4} = {01, 11, 11}. Then we have v1=0, v2=1, v3=0, and v4=2. To clarify the notation ςi(r), we now explain how to calculate ς1(1). We arbitrarily number bits from left to right; that is, in any given bit string, bit 1 is the leftmost bit and bit 2 is the rightmost bit. From the definition of ςi(r) we see that
(5)ς1(1)={j:xj(1)=x1(1)}.
Since x1 = 00, we see that x1(1) = 0 (i.e., the leftmost bit is 0). Then (5) can be written as
(6)ς1(1)={j:xj(1)=0}.
But xj(1) = 0 for xj∈{00,01}, which in turn indicates that xj(1) = x1(1) for j∈[1,2]; therefore, ς1(1) = {1, 2}. Continuing this process, we see that
(7)ς1(1)={1,2},ς1(2)={1,3},ς2(1)={1,2},ς2(2)={2,4},ς3(1)={3,4},ς3(2)={1,3},ς4(1)={3,4},ς4(2)={2,4}.Mutation. Mutation operates independently on each candidate solution by probabilistically reversing each bit in each candidate solution. Suppose that the event that each bit of a candidate solution is flipped is stochastically independent and occurs with probability pm∈(0,1). Then the probability that candidate solution xi mutates to become xj can be written as
(8)Pr(xi⟶xj)=pmHij(1-pm)q-Hij>0,
where q is the number of bits in each candidate solution and Hij represents the Hamming distance between bit strings xi and xj.


## 2.1. Overview of Biogeography-Based Optimization

This section provides an overview of BBO. The review in this section is very general because it applies to optimization problems with either real domains, integer domains, binary domains, or combinations thereof.BBO is a new optimization approach, inspired by biogeography theory, to solve general optimization problems. A biogeography habitat corresponds to acandidate solution to the optimization problem. A multiset of biogeography habitats corresponds to apopulation of candidate solutions. Habitat suitability index (HSI) in biogeography corresponds to the goodness of a candidate solution, which is also calledfitnessin standard EA notation. Like other EAs [27], BBO probabilistically shares information between candidate solutions to improve candidate solution fitness. In BBO, each candidate solution is comprised of a set of features, which are also called independent variables or decision variables in the optimization literature. Note that the decision variables can be taken from sets of real numbers, integers, binary numbers, or combinations thereof, depending on the problem. Each candidate solution immigrates features from other candidate solutions based on its immigration rate and emigrates features to other candidate solutions based on its emigration rate. BBO consists of two main steps: migration and mutation.Migration. Migration is a probabilistic operator that is intended to improve a candidate solution yk. For each feature of a given candidate solution yk, the candidate solution’s immigration rate λk is used to probabilistically decide whether or not to immigrate. If immigration occurs, then the emigrating candidate solution yj is probabilistically chosen based on the emigration rate μj. Migration is written as follows:
(1)yk(r)⟵yj(r),
where r is a candidate solution feature. In BBO, each candidate solution yk has its own immigration rate λk and emigration rate μk. A good candidate solution has relatively high μ and low λ, while the converse is true for a poor candidate solution. According to [1], the immigration rate λk and emigration rate μk of the candidate solution yk are based on linear functions and are calculated as
(2)λk=1-fitness(yk)μk=fitness(yk),
where fitness denotes candidate solution fitness value, which is normalized to the range [0, 1]. The probabilities of immigrating to yk and of emigrating from yk are calculated as
(3)Pr(immigrationtoyk)=λkPr(emigrationfromyk)=μk∑j=1Nμj,
where N is the population size.Mutation. Mutation is a probabilistic operator that randomly modifies a candidate solution feature. The purpose of mutation is to increase diversity among the population, just as in other EAs. Mutation of the ith candidate solution is implemented as shown in Algorithm 1.Algorithm 1
F
or each candidate solution feature
r  inyi
Ifrand(0,1)  <pm
yi(r)←  rand(Lr,Ur)
End if
Next
rIn the mutation logic (Algorithm1), rand (a,b) is a uniformly distributed random number between a and b, pm is the mutation rate, and Lr and Ur are the lower and upper search bounds of the rth independent variable. The above logic mutates each independent variable with a probability of pm. If mutation occurs for a given independent variable, then that independent variable is replaced with a random number within its search domain.A description of one generation of BBO is given in Algorithm2.Algorithm 2:One generation of the BBO algorithm. y is the entire population of candidate solutions, yk is the kth candidate solution, and yk(r) is the rth feature of yk.
For eachyk, define emigration rate μk proportional to fitness of yk, μk∈[0,1]
For eachyk, define immigration rate λk=1-μk
For eachyi
For each candidate solution featurer
Use λi to probabilistically decide whether to immigrate to yi (see (2) and (3))
If immigrating then
Use{μk} to probabilistically select the emigrating solution yk (see (2) and (3))
yi(r)←yk(r)
End if
Next candidate solution feature
Probabilistically decide whether to mutateyi
Next solution


## 2.2. Transition Probability of Biogeography-Based Optimization

Two main steps in BBO are significant: migration and mutation. The transition probability of one generation includes migration probability and mutation probability. Consider a problem whose search space{x1,…,xn} is binary. The search space is the set of all bit strings xi consisting of q bits each. Therefore, the cardinality of the search space is n=2q. Suppose that BBO is currently in the tth generation. Based on the previously derived Markov chain model for BBO [11], the probability that migration results in the kth candidate solution yk are equal to xi at generation t+1 is given by
(4)Pr(yk,t+1=xi)=∏r=1q[(1-λk)10(yk(r)-xi(r))+λk∑j∈ςi(r)vjμj∑j=1nvjμj],
where10 is the indicator function on the set 0 (i.e.,10(a)=1 if a = 0, and10(a) = 0 if a≠0), r denotes the index of the candidate solution feature (i.e., the bit number), λk denotes the immigration rate of candidate solution yk, μj denotes the emigration rate of candidate solution xj, and vj denotes the number of xj individuals in the population. The notation ςi(r) in (4) denotes the set of search space indices j such that the rth bit of xj is equal to the rth bit of xi. That is, ςi(r) = {j:xj(r)=xi(r)}. Note that the first term in the product on the right side of (4) denotes the probability that yk,t+1(r)=xi(r) if immigration of the rth candidate solution feature does not occur, and the second term denotes the probability if immigration of the rth candidate solution feature does occur.Example 1.
To clarify the notations in (4), an example is presented. We use the notation vj to denote the number of xj individuals in the population. Suppose we have a two-bit problem (q=2, n=4) with a population size N=3. The search space consists of the bit strings x={x1, x2, x3, x4} = {00, 01, 10, 11}. Suppose that the individuals in the current population are y={x2, x4, x4} = {01, 11, 11}. Then we have v1=0, v2=1, v3=0, and v4=2. To clarify the notation ςi(r), we now explain how to calculate ς1(1). We arbitrarily number bits from left to right; that is, in any given bit string, bit 1 is the leftmost bit and bit 2 is the rightmost bit. From the definition of ςi(r) we see that
(5)ς1(1)={j:xj(1)=x1(1)}.
Since x1 = 00, we see that x1(1) = 0 (i.e., the leftmost bit is 0). Then (5) can be written as
(6)ς1(1)={j:xj(1)=0}.
But xj(1) = 0 for xj∈{00,01}, which in turn indicates that xj(1) = x1(1) for j∈[1,2]; therefore, ς1(1) = {1, 2}. Continuing this process, we see that
(7)ς1(1)={1,2},ς1(2)={1,3},ς2(1)={1,2},ς2(2)={2,4},ς3(1)={3,4},ς3(2)={1,3},ς4(1)={3,4},ς4(2)={2,4}.Mutation. Mutation operates independently on each candidate solution by probabilistically reversing each bit in each candidate solution. Suppose that the event that each bit of a candidate solution is flipped is stochastically independent and occurs with probability pm∈(0,1). Then the probability that candidate solution xi mutates to become xj can be written as
(8)Pr(xi⟶xj)=pmHij(1-pm)q-Hij>0,
where q is the number of bits in each candidate solution and Hij represents the Hamming distance between bit strings xi and xj.


## 3. Convergence of Biogeography-Based Optimization

The previous section reviewed the BBO algorithm and its Markov model. In this section, which comprises the main contribution of this paper, we use the results of Section2 to analyze the convergence behavior of BBO. Section 3.1 gives some basic foundations of Markov transition matrices, including notation and basic theorems that we will need later. Section 3.2 reviews previously published Markov theory as it relates to BBO. This leads to Section 3.3, which obtains some important properties of the BBO transition matrix. Section 3.4 uses transition matrix properties to analyze BBO convergence to the solution of a global optimization problem. This leads to Section 3.5, which uses the BBO convergence analysis to obtain an estimate of the convergence rate.


### 3.1. Preliminary Foundations of Markov Theory

A finite Markov chain is a random process which has a finite number of possible state valuesS={si} (i=1,2,…,T), where T is the total number of states, which is the cardinality |S|. The probability that the system transitions from state si to sj at time step t is given by pij(t), which is called the transition probability. The T×T matrix P = (pij (t)) is called the transition matrix, where pij∈[0,1] for i, j∈[1,T], and ∑j=1Tpij=1 for all i. The P matrix is calledstochastic because the elements in each row sum to 1. If the transition probability is independent of t, that is, pij(t1)=pij(t2) for all i,  j∈ [1, T] and for all t1 and t2, then the Markov chain is said to behomogeneous. Given an initial probability distribution of states π(0) as a row vector, the probability distribution of the Markov chain after t steps is given by π(t)=π(0)Pt. Therefore, a homogeneous finite Markov chain is completely specified by π(0) and P, and the limiting distribution as t→∞ depends on the structure of P. For homogeneous finite Markov chains, we have the following two theorems [28, 29].Theorem 2 (see [28, page 123]).
LetP be a primitive stochastic matrix of order T; that is, all of the elements of Pk are positive for some integer k. Then Pk converges as k→∞ to a stochastic matrix which has all nonzero entries. That is, for all i, j∈[1,T],
(9)P∞=lim⁡k→∞(pij)(k)=(π⋮π)T×T,
where π=(π1,…,πT), and πj≠0 for 1≤j≤T.Theorem 3 (see [28, page 126]).
LetP be a stochastic matrix of order T with the structure P=[C0RQ], where C is a primitive stochastic matrix of order m and R, Q≠0. Then Pk converges as k→∞ to a stochastic matrix. That is,
(10)P∞=lim⁡k→∞(pij)(k)=(π⋮π)T×T,
where π=(π1,…,πm,0,…,0), and πj≠0 for 1≤j≤m<T.We will use these theorems in Section3.3 to derive important properties of the BBO transition matrix and in Section 3.4 to derive BBO convergence properties.


### 3.2. BBO Markov Theory

In previous work [11], the transition probability of BBO with migration and mutation was obtained. This provides us with the probability Pr(u∣v) of transitioning in one generation from population vector v = [v1,v2,…,vn], where vi is the number of candidate solutions xi in the population and n is the size of the search space, to population vector u=[u1,u2,…,un]. BBO can be described as a homogeneous finite Markov chain: the state of BBO is defined as the population vector, so the element pij of the state transition matrix P is obtained by computing Pr(u∣v) for each possible v and each possible u. Namely, pij=Pr(u∣v) denotes the probability that the ith population vector, denoted as v, transitions to the jth population vector, denoted as u, where i,  j∈ [1, T]. Note that the cardinality of the state space is |S|=T, where T is the total number of possible populations; that is, T is the number of possible vectors v and the number of possible vectors u. The number T can be calculated in several different ways, as discussed in [11].Let theN×n matrix  M=(mkj) and the n×n matrix  U=(uij) be intermediate transition matrices corresponding to only migration and only mutation, respectively, where N is the population size and n is the cardinality of the search space. Note that mkj=Pr(yk,t+1=xj)≥0 and uij=Pr(xj→xi)>0. That is, mkj is the probability that the kth individual in the population transitions to the jth individual in the search space when only migration is considered, and uij is the probability that the jth individual in the search space transitions to the ith individual in the search space when only mutation is considered. We can use [11] to obtain the transition probability of the ith population state vector v to the jth population state vector u as
(11)pij=Pr(u∣v)=∑J∈Y∏k=1N∏i=1n[Qki(v)]Jki,whereY≡{∑KNJ∈RN×n:Jki∈{0,1},mmmmmmm∑i=1nJki=1∀k,∑k=1NJki=ui∀i},
where Qki(v) is a single element of the product of M and U. The matrix composed of the  Qki(v) elements can be represented as [Qki(v)] = MUT, where k∈[1,N] and i∈[1,n]. Qki(v) denotes the probability that the kth migration trial followed by mutation results in candidate solution xi. Note that Qki(v) is a scalar, and transition matrix P = (pij) is a T×T matrix, each element of which can be obtained by (11).Example 4.
Here we use a simple example based on [11] to clarify (11). Consider a simple BBO experiment in which a trial of migration and mutation can result in one of four possible outcomes x1, x2, x3, and x4 with probabilities Qi1, Qi2, Qi3, and Qi4, respectively. Index i refers to the migration trial number (i.e., the “For each yi” loop in Algorithm 2). Assume that the total number of trials (i.e., the population size N) is equal to 2. Suppose that the probabilities are given as follows:Trial 1:Q11=0.1, Q12=0.3, Q13=0.5, Q14=0.1;Trial 2:Q21=0.1, Q22=0.1, Q23=0.1, Q24=0.7.In this example, we calculate the probability that x1 and x4 occur after two migration trials. In order to calculate this probability, let c=[C1,C2,C3,C4] denote a vector of random variables, where Ci is the total number of times that xi occurs after two migration trials. Based on (11), we use C1=1, C2=0, C3=0, and C4=1 to obtain
(12)Pr(C1=1,C2=0,C3=0,C4=1)=∑J∈Y∏k=12∏i=14[Qki]Jki,
where
(13)Y={∑i=14Jki=1∀k,∑k=12Jki=Ci∀iJ∈R2×4:Jki∈{0,1},mmm∑i=14Jki=1∀k,∑k=12Jki=Ci∀i}.
According to (13), J belongs to Y if it satisfies the following conditions: (1) J is a 2×4 matrix; (2) each element of  J is either 0 or 1; (3) the elements in each row of  J add up to 1; and (4) the elements in the ith column of J add up to Ci.Note from [11] that the cardinality of Y is
(14)|Y|=N!(C1!⋯Cn!).
The number of matrices J(t) that satisfy these conditions is calculated as N!/(C1!⋯Cn!)=2!/(1!1!0!1!)=2, and the J(t) matrices are found as follows:
(15)J(1)=[10000001],J(2)=[00011000].
Substituting these matrices into (12) gives
(16)Pr(C1=1,C2=0,C3=0,C4=1)=∑t=12(Q11J11(t)Q12J12(t)Q13J13(t)Q14J14(t))×(Q21J21(t)Q22J22(t)Q23J23(t)Q24J24(t))=Q11Q24+Q14Q21=0.08.


### 3.3. BBO Transition Matrix Properties

Recall that the migration probability and mutation probability can be calculated by (3) and (8), respectively: mkj = Pr (yk,t+1=xj) ≥ 0 and uij = Pr (xj→xi) > 0. Therefore, M is a nonnegative stochastic matrix; although it is not a transition matrix since it is not square, each row sums to 1. We also see that U is a positive left stochastic matrix; that is, each of its columns sums to 1. We now present two theorems that show that there is a nonzero probability of obtaining any individual in the search space from any individual in a BBO population after migration and mutation. This means that there is a nonzero probability of transitioning from any population vector u to any other population vector v in one generation, which means that the BBO transition matrix is primitive.Theorem 5.
IfM is a positive stochastic matrix and U is a positive left stochastic matrix, then the product MUT is positive.Proof.
IfM is positive and stochastic, then every entry of M is positive; that is, mkj > 0 for k∈[1,N] and j∈[1,n], and ∑j=1nmkj=1 for all k. Similarly, if U is positive, then every entry of U is positive; that is, uij>0 for i,j∈[1,n]. Therefore, by matrix multiplication, (MUT)ki=∑j=1nmkjuij>0 for k∈[1,N] and i∈[1,n].Theorem 6.
The transition matrix of BBO with migration and mutation is primitive.Proof.
From (11) we know that if Qki(v) = [MUT]ki>0 for all k∈[1,N] and i∈[1,n], then pij=∑J∈Y∏k=1N∏i=1n[Qki(v)]Jki>0 for i,j∈[1,T], where Y is given in (11). So the transition matrix P=(pij) of BBO is positive. Therefore, P is primitive since every positive transition matrix is primitive.Corollary 7.
There exists a unique limiting distribution for the states of the BBO Markov chain. Also, the probability that the Markov chain is in theith state at any time is nonzero for all i∈[1,T].Proof.
Corollary7 is an immediate consequence of Theorems 2 and 6.


### 3.4. Convergence Properties of BBO

Before we obtain the convergence properties of BBO, some precise definitions of the termconvergence are required [15]. Assume that the search space of a global optimization problem is I with cardinality |I|=n. Further assume that the BBO algorithm with population size N consists of both migration and mutation, as shown in Algorithm 2.Definition 8.
LetA(t)={ai(t)∣i∈[1,N],ai(t)∈I} be the population at generation t, where N is the population size, ai(t) denotes an individual representing a candidate solution in the search space I, and A(t) may contain duplicate elements; f:I→R denotes a fitness function assigning real values to individuals; I*={a*|a*=argmax⁡{f(a)∣a∈I}} is a subset in the search space, each member of which has the globally maximum fitness; and the best individuals in the population at generation t are I*(t)={aj*(t)}⊂A(t), where f(aj*(t))≥f(ai(t)) for all aj*(t)∈I*(t) and for all i∈[1,N].We use the notationa*(t) to denote an arbitrary element of I*(t) (i.e., one of the best individuals in the population at generation t). Because of migration and mutation, a*(t) and its fitness will change randomly over time. As t→∞, the convergence, or lack of convergence, of a*(t) to the subset I* indicates whether or not the BBO algorithm is globally convergent. That is, BBO is said to converge if
(17)Pr(lim⁡t→∞a*(t)∈I*)=1⟺Pr(a*∈lim⁡t→∞A(t))=1.
Note that a*(t) is not necessarily unique. However, Definition 8 states that the BBO algorithm is globally convergent if and only if (17) holds for every a*(t). Clearly, the evolution of a*(t) is a homogeneous finite Markov chain, which we call an a*(t)-chain.Now we sort all the states ofI in order of descending fitness; that is, I={I1,…,In}, and f(I1)≥f(I2)≥⋯≥f(In). We define S as the set of indices of I; that is, S={1,2,…,n}. Further, we define S* as the elements {j} of S such that Ij∈I*; that is, Ij∈I* for all j∈S*. This leads to the following definition.Definition 9.
LetP^=(p^ij) be the transition matrix of an a*(t)-chain, where p^ij for i,j∈[1,n] is the probability that a*(t)=Ii transitions to a*(t+1)=Ij. The BBO algorithm converges to a global optimum if and only if a*(t) transitions from any state i∈I to I* as t→∞ with probability one, that is, if
(18)lim⁡t→∞∑j∈S*(P^t)ij=1,∀i∈S.
As noted earlier, there may be more than one a*(t)-chain since more than one element of the search space may have a globally maximum fitness. Definition 9 states that the BBO algorithm converges to a global optimum if and only if (18) holds for every a*(t)-chain. Also note that P^ depends on the other individuals in the population at generation t. Definition 9 states that the BBO algorithm converges to a global optimum if and only if (18) holds for every possible P^ transition matrix for every a*(t)-chain.Theorem 10.
If the transition matrixP^=(p^ij) of an a*(t)-chain is a positive stochastic matrix, then BBO with migration and mutation does not converge to any of the global optima.Proof.
Since every positive matrix is also a primitive one, it follows by Theorem2 that the limiting distribution of P^ is unique with all nonzero entries. Therefore, for any i∈S,
(19)lim⁡t→∞∑j∈S*(P^t)ij=1-lim⁡t→∞∑j∈S-S*(P^t)ij=1-∑j∈S-S*lim⁡t→∞(P^t)ij=1-∑j=|S*|+1|S|πj<1,
where we use the notation S-S* to denote all elements of S that do not belong to S*. We see that (18) is not satisfied, which completes the proof.Theorem 11.
If the transition matrixP^=(p^ij) of an a*(t)-chain is a stochastic matrix with the structure P^=[C0RQ], where C is a positive stochastic matrix of order |S*|, and R, Q≠0, then the BBO algorithm converges to one or more of the global optima.Proof.
From Theorem3, we see that, for all i, j∈S,
(20)P^∞=lim⁡t→∞(p^ij)t=(π⋮π)|S|×|S|,
where π=(π1,…,π|S*|,0,…,0), πj≠0 for 1≤j≤|S*|, and ∑j=1|S*|πj=1. It follows directly that, for any i∈S,
(21)lim⁡t→∞∑j∈S*(P^t)ij=∑j∈S*lim⁡t→∞(P^t)ij=∑j=1|S*|πj=1.
We see that (18) is satisfied which completes the proof.Theorems10 and 11 can be applied directly to determine the global convergence of BBO if the structure of the transition matrix of the Markov chain can be determined, as we will show in the remainder of this section. In particular, we will formalize the observation that the transition matrix of BBO without elitism satisfies the conditions of Theorem 10 (as stated in Theorem 6). We will further show that the transition matrix of the a*(t)-chain of BBO with elitism satisfies the conditions of Theorem 11.Elitism. We now discuss a modified BBO which useselitism, an idea which is also implemented in many EAs. There are many ways to implement elitism, but here we define elitism as the preservation of the best individual at each generation in a separate partition of the population space. This enlarges the population size by one individual; the elite individual increases the population size from N to N+1. However, note that the population size is still constant (i.e., equal to N+1) from one generation to the next. The elite individual does not take part in recombination or mutation but is maintained separately from the other N-members of the population. At each generation, if an individual in the N-member main population is better than the elite individual, then the elite individual is replaced with a copy of the better individual.Relative to a standardN-member BBO population, elite BBO increases the number of possible population distributions by a factor of n, which is the search space size. That is, each possible population distribution of the N-member main population could also include one of n elite individuals. The number of possible population distributions increases by a factor of n, from T (see Section 3.1) to nT. We order these new states so that each group of n states has the same elite individual. Also, the elite individual in the mth group of n states is the mth best individual in the search space, for m=1,…,n.The elitist-preserving process can be represented by an upgrade transition matrixO, which contains the probabilities that each population distribution of the (N+1)-member population transitions to some other population distribution after the elitist-preserving step. That is, the element in the ith row and jth column of O, denoted as O(i,j), is the probability that the ith population distribution transitions to the jth population distribution after the step in which the elite individual is replaced with the best individual from the N-member main population. The upgrade matrix is similar to the one in [29]; it does not include the effects of migration or mutation but only includes the elitism-preserving step. The upgrade matrix only includes the probability of changing the elite individual; it does not include the probability of changing the N-member main population, since it does not include migration or mutation. If there are no individuals in the N-member main population that are better than the elite individual, then the elite individual does not change. The structure of the upgrade matrix O can be written as
(22)O=(O110⋯0O21O22⋱⋮⋮⋱⋱0On1⋯On,n-1Onn),
where each Oij matrix is T×T, where T is the number of population distributions in an EA with a population size of N and search space cardinality of n. O11 is the identity matrix since the first n population distributions have the global optimum as their elite individual, and the elite individual can never be improved from the global optimum. Matrices Oaa with a≥2 are diagonal matrices composed of all zeros and ones. Since the population distributions are ordered by grouping common elite individuals and since elite individuals in the population distribution ordering are in order of decreasing fitness, the super block diagonals in O are zero matrices as shown in (22); that is, there is zero probability that the ith population distribution transitions to the jth population distribution if i<j. So the Markov chain of elite BBO can be described as
(23)P+=diag⁡(P)O=(P0⋯00P⋱⋮⋮⋱⋱00⋯0P)(O110⋯0O21O22⋱⋮⋮⋱⋱0On1⋯On,n-1Onn)=(PO110⋯0PO21PO22⋱⋮⋮⋱⋱0POn1⋯POn,n-1POnn),
where P is the T×T transition matrix described in Section 3.1.Example 12.
To explain the update matrixO described in (22), a simple example is presented. Suppose there exists a search space consisting of n=3 individuals which are x={x1,x2,x3} where the fitness of x1 is the lowest and the fitness of x3 is the highest. Suppose the main population size is N=1, so the elitist population size is N+1=2. So there are nine possible populations before the elitist-preserving step, and they are {C1,C2,C3,C4,C5,C6,C7,C8,C9} = {{x3,x3}, {x3,x2}, {x3,x1}, {x2,x3}, {x2,x2}, {x2,x1}, {x1,x3}, {x1,x2}, {x1,x1}}. Note that the first element in each population is the elite individual and the last N element (N=1 in this example) is the main population. Also note that the populations are ordered in such a way that the first three have the most fit individual as their elite individual, the next three have the second most fit individual as their elite individual, and the last three have the least fit individual as their elite individual. The update matrix O is a 9×9 matrix.The populationC1={x3,x3} transitions to the population C1={x3,x3} with probability 1; that is, O(1,1) = 1. Population C1 cannot transition to any other population Ci(i≠1); that is, O(1,i) = 0 for i≠1. Similarly, population C2={x3,x2} transitions to C2 with probability 1 since the elite x3 is better than the main-member population x2; therefore, O(2,2)=1, and O(2,i)=0 for i≠2. Continuing with this reasoning, we obtain the O matrix as follows:
(24)O=[O11O12O13O21O22O23O31O32O33]=[101001100000001000001100000000001000000000001],
where each Oij matrix is 3×3 and the blank elements are 0.Now we consider the convergence of thea*(t)-chain, which is the sequence of elite individuals in the elite BBO algorithm. If the elite individual is equal to the global optimum, we call this an absorbing state of the a*(t)-chain. Recall that the elite individual in elite BBO can only be replaced by one with better fitness. Therefore, the a*(t)-chain of elite BBO contains three classes of states: (1) at least one absorbing state, (2) nonabsorbing states which transition to absorbing states in one step, and (3) nonabsorbing states which transition to nonabsorbing states in one step. So the transition matrix P^ of the a*(t)-chain, which we introduced in (18)–(21), can be written as
(25)P^=(Ik0RQ),
where Ik is a k×k unit matrix corresponding to optimal individuals (k is the number of optima), R is a matrix of order (|S|-k)×k corresponding to nonabsorbing states that transition to absorbing states (|S| is the cardinality of the state space S, so |S|-k is the number of nonabsorbing states), and Q is a matrix of order (|S|-k)×(|S|-k) corresponding to nonabsorbing states that transition to nonabsorbing states. The matrix P^ of (25) has the same structure as the matrix P^ in Theorem 11. It follows from Theorem 11 that the a*(t)-chain of elite BBO is globally convergent.These results are similar to the canonical GA [29], which is proven to never converge to the global optimum, but elitist variants of which are proven to converge to the global optimum. We sum up these results in the following corollary.Corollary 13.
BBO with migration and mutation does not converge to any of the global optima, but elite BBO, which preserves the best individual at each generation, converges to the global optimum.Proof.
This is an immediate consequence of Theorems6 and 10 (the nonconvergence of BBO without elitism), Theorem 11 (the convergence of BBO with elitism), and the discussion above.


### 3.5. Convergence Rate

The previous subsection analyzed the convergence properties of elite BBO, and this subsection discusses its convergence rate. The transition matrix of elite BBO aftert steps can be found from (25) as follows:
(26)P^t=(Ik0NtRQt),
where Nt=I+Q+Q2+⋯+Qt-1=(I-Q)-1(I-Qt). If ∥Q∥<1 the limiting distribution of the Markov chain of BBO can be found from P^∞, which can be written as
(27)P^∞=lim⁡t→∞P^t=[Ik0(I-Q)-1R0].
Modified BBO with elitism has been proven to converge to a global optimum in the previous subsection, and there exists a limiting distribution π*=π(0)P^∞, where π(0) = [π1(0), …, πk(0), πk+1(0),…, π|s|(0)]=[π-1(0),π-2(0)] (recall that k is the number of global optima). The convergence rate estimate of elite BBO can be obtained as follows.Theorem 14.
If∥Q∥=ρ<1, the convergence rate of elite BBO satisfies ∥π(t)-π*∥≤O(ρt).Proof.
Consider(28)∥π(t)-π*∥=∥π(0)P^t-π(0)P^∞∥=∥π(0)[[Ik0NtRQt]-[Ik0(I-Q)-1R0]]∥=∥[π-1(0)π-2(0)][00(Nt-(I-Q)-1)RQt]∥=∥π-2(0)[(Nt-(I-Q)-1)RQt]∥≤∥π-2(0)∥(∥Nt-(I-Q)-1∥·∥R∥+∥Qt∥)=∥π-2(0)∥(∥(I-Q)-1(I-Qt)-(I-Q)-1∥mmmmmmm·∥R∥+∥Qt∥)=∥π-2(0)∥(∥(I-Q)-1-(I-Q)-1Qt-(I-Q)-1∥mmmmmmm·∥R∥+∥Qt∥)≤∥π-2(0)∥(∥-(I-Q)-1∥·∥R∥+1)∥Qt∥≤∥π-2(0)∥(∥R∥1-∥Q∥+1)∥Q∥t=O(ρt).
Note that elite BBO is guaranteed to converge to a global optimum regardless of the initial state. In addition, note that we can improve the convergence rate bound by decreasing the parameter ρ. That is, reducing the number of nonabsorbing states which transition to other nonabsorbing states can accelerate the convergence of elite BBO. In spite of differences between GAs and BBO [12], we see from Theorem 14 that the convergence rate of BBO with elitism is very similar to that of GAs [22, Theorem 5].


## 3.1. Preliminary Foundations of Markov Theory

A finite Markov chain is a random process which has a finite number of possible state valuesS={si} (i=1,2,…,T), where T is the total number of states, which is the cardinality |S|. The probability that the system transitions from state si to sj at time step t is given by pij(t), which is called the transition probability. The T×T matrix P = (pij (t)) is called the transition matrix, where pij∈[0,1] for i, j∈[1,T], and ∑j=1Tpij=1 for all i. The P matrix is calledstochastic because the elements in each row sum to 1. If the transition probability is independent of t, that is, pij(t1)=pij(t2) for all i,  j∈ [1, T] and for all t1 and t2, then the Markov chain is said to behomogeneous. Given an initial probability distribution of states π(0) as a row vector, the probability distribution of the Markov chain after t steps is given by π(t)=π(0)Pt. Therefore, a homogeneous finite Markov chain is completely specified by π(0) and P, and the limiting distribution as t→∞ depends on the structure of P. For homogeneous finite Markov chains, we have the following two theorems [28, 29].Theorem 2 (see [28, page 123]).
LetP be a primitive stochastic matrix of order T; that is, all of the elements of Pk are positive for some integer k. Then Pk converges as k→∞ to a stochastic matrix which has all nonzero entries. That is, for all i, j∈[1,T],
(9)P∞=lim⁡k→∞(pij)(k)=(π⋮π)T×T,
where π=(π1,…,πT), and πj≠0 for 1≤j≤T.Theorem 3 (see [28, page 126]).
LetP be a stochastic matrix of order T with the structure P=[C0RQ], where C is a primitive stochastic matrix of order m and R, Q≠0. Then Pk converges as k→∞ to a stochastic matrix. That is,
(10)P∞=lim⁡k→∞(pij)(k)=(π⋮π)T×T,
where π=(π1,…,πm,0,…,0), and πj≠0 for 1≤j≤m<T.We will use these theorems in Section3.3 to derive important properties of the BBO transition matrix and in Section 3.4 to derive BBO convergence properties.


## 3.2. BBO Markov Theory

In previous work [11], the transition probability of BBO with migration and mutation was obtained. This provides us with the probability Pr(u∣v) of transitioning in one generation from population vector v = [v1,v2,…,vn], where vi is the number of candidate solutions xi in the population and n is the size of the search space, to population vector u=[u1,u2,…,un]. BBO can be described as a homogeneous finite Markov chain: the state of BBO is defined as the population vector, so the element pij of the state transition matrix P is obtained by computing Pr(u∣v) for each possible v and each possible u. Namely, pij=Pr(u∣v) denotes the probability that the ith population vector, denoted as v, transitions to the jth population vector, denoted as u, where i,  j∈ [1, T]. Note that the cardinality of the state space is |S|=T, where T is the total number of possible populations; that is, T is the number of possible vectors v and the number of possible vectors u. The number T can be calculated in several different ways, as discussed in [11].Let theN×n matrix  M=(mkj) and the n×n matrix  U=(uij) be intermediate transition matrices corresponding to only migration and only mutation, respectively, where N is the population size and n is the cardinality of the search space. Note that mkj=Pr(yk,t+1=xj)≥0 and uij=Pr(xj→xi)>0. That is, mkj is the probability that the kth individual in the population transitions to the jth individual in the search space when only migration is considered, and uij is the probability that the jth individual in the search space transitions to the ith individual in the search space when only mutation is considered. We can use [11] to obtain the transition probability of the ith population state vector v to the jth population state vector u as
(11)pij=Pr(u∣v)=∑J∈Y∏k=1N∏i=1n[Qki(v)]Jki,whereY≡{∑KNJ∈RN×n:Jki∈{0,1},mmmmmmm∑i=1nJki=1∀k,∑k=1NJki=ui∀i},
where Qki(v) is a single element of the product of M and U. The matrix composed of the  Qki(v) elements can be represented as [Qki(v)] = MUT, where k∈[1,N] and i∈[1,n]. Qki(v) denotes the probability that the kth migration trial followed by mutation results in candidate solution xi. Note that Qki(v) is a scalar, and transition matrix P = (pij) is a T×T matrix, each element of which can be obtained by (11).Example 4.
Here we use a simple example based on [11] to clarify (11). Consider a simple BBO experiment in which a trial of migration and mutation can result in one of four possible outcomes x1, x2, x3, and x4 with probabilities Qi1, Qi2, Qi3, and Qi4, respectively. Index i refers to the migration trial number (i.e., the “For each yi” loop in Algorithm 2). Assume that the total number of trials (i.e., the population size N) is equal to 2. Suppose that the probabilities are given as follows:Trial 1:Q11=0.1, Q12=0.3, Q13=0.5, Q14=0.1;Trial 2:Q21=0.1, Q22=0.1, Q23=0.1, Q24=0.7.In this example, we calculate the probability that x1 and x4 occur after two migration trials. In order to calculate this probability, let c=[C1,C2,C3,C4] denote a vector of random variables, where Ci is the total number of times that xi occurs after two migration trials. Based on (11), we use C1=1, C2=0, C3=0, and C4=1 to obtain
(12)Pr(C1=1,C2=0,C3=0,C4=1)=∑J∈Y∏k=12∏i=14[Qki]Jki,
where
(13)Y={∑i=14Jki=1∀k,∑k=12Jki=Ci∀iJ∈R2×4:Jki∈{0,1},mmm∑i=14Jki=1∀k,∑k=12Jki=Ci∀i}.
According to (13), J belongs to Y if it satisfies the following conditions: (1) J is a 2×4 matrix; (2) each element of  J is either 0 or 1; (3) the elements in each row of  J add up to 1; and (4) the elements in the ith column of J add up to Ci.Note from [11] that the cardinality of Y is
(14)|Y|=N!(C1!⋯Cn!).
The number of matrices J(t) that satisfy these conditions is calculated as N!/(C1!⋯Cn!)=2!/(1!1!0!1!)=2, and the J(t) matrices are found as follows:
(15)J(1)=[10000001],J(2)=[00011000].
Substituting these matrices into (12) gives
(16)Pr(C1=1,C2=0,C3=0,C4=1)=∑t=12(Q11J11(t)Q12J12(t)Q13J13(t)Q14J14(t))×(Q21J21(t)Q22J22(t)Q23J23(t)Q24J24(t))=Q11Q24+Q14Q21=0.08.


## 3.3. BBO Transition Matrix Properties

Recall that the migration probability and mutation probability can be calculated by (3) and (8), respectively: mkj = Pr (yk,t+1=xj) ≥ 0 and uij = Pr (xj→xi) > 0. Therefore, M is a nonnegative stochastic matrix; although it is not a transition matrix since it is not square, each row sums to 1. We also see that U is a positive left stochastic matrix; that is, each of its columns sums to 1. We now present two theorems that show that there is a nonzero probability of obtaining any individual in the search space from any individual in a BBO population after migration and mutation. This means that there is a nonzero probability of transitioning from any population vector u to any other population vector v in one generation, which means that the BBO transition matrix is primitive.Theorem 5.
IfM is a positive stochastic matrix and U is a positive left stochastic matrix, then the product MUT is positive.Proof.
IfM is positive and stochastic, then every entry of M is positive; that is, mkj > 0 for k∈[1,N] and j∈[1,n], and ∑j=1nmkj=1 for all k. Similarly, if U is positive, then every entry of U is positive; that is, uij>0 for i,j∈[1,n]. Therefore, by matrix multiplication, (MUT)ki=∑j=1nmkjuij>0 for k∈[1,N] and i∈[1,n].Theorem 6.
The transition matrix of BBO with migration and mutation is primitive.Proof.
From (11) we know that if Qki(v) = [MUT]ki>0 for all k∈[1,N] and i∈[1,n], then pij=∑J∈Y∏k=1N∏i=1n[Qki(v)]Jki>0 for i,j∈[1,T], where Y is given in (11). So the transition matrix P=(pij) of BBO is positive. Therefore, P is primitive since every positive transition matrix is primitive.Corollary 7.
There exists a unique limiting distribution for the states of the BBO Markov chain. Also, the probability that the Markov chain is in theith state at any time is nonzero for all i∈[1,T].Proof.
Corollary7 is an immediate consequence of Theorems 2 and 6.


## 3.4. Convergence Properties of BBO

Before we obtain the convergence properties of BBO, some precise definitions of the termconvergence are required [15]. Assume that the search space of a global optimization problem is I with cardinality |I|=n. Further assume that the BBO algorithm with population size N consists of both migration and mutation, as shown in Algorithm 2.Definition 8.
LetA(t)={ai(t)∣i∈[1,N],ai(t)∈I} be the population at generation t, where N is the population size, ai(t) denotes an individual representing a candidate solution in the search space I, and A(t) may contain duplicate elements; f:I→R denotes a fitness function assigning real values to individuals; I*={a*|a*=argmax⁡{f(a)∣a∈I}} is a subset in the search space, each member of which has the globally maximum fitness; and the best individuals in the population at generation t are I*(t)={aj*(t)}⊂A(t), where f(aj*(t))≥f(ai(t)) for all aj*(t)∈I*(t) and for all i∈[1,N].We use the notationa*(t) to denote an arbitrary element of I*(t) (i.e., one of the best individuals in the population at generation t). Because of migration and mutation, a*(t) and its fitness will change randomly over time. As t→∞, the convergence, or lack of convergence, of a*(t) to the subset I* indicates whether or not the BBO algorithm is globally convergent. That is, BBO is said to converge if
(17)Pr(lim⁡t→∞a*(t)∈I*)=1⟺Pr(a*∈lim⁡t→∞A(t))=1.
Note that a*(t) is not necessarily unique. However, Definition 8 states that the BBO algorithm is globally convergent if and only if (17) holds for every a*(t). Clearly, the evolution of a*(t) is a homogeneous finite Markov chain, which we call an a*(t)-chain.Now we sort all the states ofI in order of descending fitness; that is, I={I1,…,In}, and f(I1)≥f(I2)≥⋯≥f(In). We define S as the set of indices of I; that is, S={1,2,…,n}. Further, we define S* as the elements {j} of S such that Ij∈I*; that is, Ij∈I* for all j∈S*. This leads to the following definition.Definition 9.
LetP^=(p^ij) be the transition matrix of an a*(t)-chain, where p^ij for i,j∈[1,n] is the probability that a*(t)=Ii transitions to a*(t+1)=Ij. The BBO algorithm converges to a global optimum if and only if a*(t) transitions from any state i∈I to I* as t→∞ with probability one, that is, if
(18)lim⁡t→∞∑j∈S*(P^t)ij=1,∀i∈S.
As noted earlier, there may be more than one a*(t)-chain since more than one element of the search space may have a globally maximum fitness. Definition 9 states that the BBO algorithm converges to a global optimum if and only if (18) holds for every a*(t)-chain. Also note that P^ depends on the other individuals in the population at generation t. Definition 9 states that the BBO algorithm converges to a global optimum if and only if (18) holds for every possible P^ transition matrix for every a*(t)-chain.Theorem 10.
If the transition matrixP^=(p^ij) of an a*(t)-chain is a positive stochastic matrix, then BBO with migration and mutation does not converge to any of the global optima.Proof.
Since every positive matrix is also a primitive one, it follows by Theorem2 that the limiting distribution of P^ is unique with all nonzero entries. Therefore, for any i∈S,
(19)lim⁡t→∞∑j∈S*(P^t)ij=1-lim⁡t→∞∑j∈S-S*(P^t)ij=1-∑j∈S-S*lim⁡t→∞(P^t)ij=1-∑j=|S*|+1|S|πj<1,
where we use the notation S-S* to denote all elements of S that do not belong to S*. We see that (18) is not satisfied, which completes the proof.Theorem 11.
If the transition matrixP^=(p^ij) of an a*(t)-chain is a stochastic matrix with the structure P^=[C0RQ], where C is a positive stochastic matrix of order |S*|, and R, Q≠0, then the BBO algorithm converges to one or more of the global optima.Proof.
From Theorem3, we see that, for all i, j∈S,
(20)P^∞=lim⁡t→∞(p^ij)t=(π⋮π)|S|×|S|,
where π=(π1,…,π|S*|,0,…,0), πj≠0 for 1≤j≤|S*|, and ∑j=1|S*|πj=1. It follows directly that, for any i∈S,
(21)lim⁡t→∞∑j∈S*(P^t)ij=∑j∈S*lim⁡t→∞(P^t)ij=∑j=1|S*|πj=1.
We see that (18) is satisfied which completes the proof.Theorems10 and 11 can be applied directly to determine the global convergence of BBO if the structure of the transition matrix of the Markov chain can be determined, as we will show in the remainder of this section. In particular, we will formalize the observation that the transition matrix of BBO without elitism satisfies the conditions of Theorem 10 (as stated in Theorem 6). We will further show that the transition matrix of the a*(t)-chain of BBO with elitism satisfies the conditions of Theorem 11.Elitism. We now discuss a modified BBO which useselitism, an idea which is also implemented in many EAs. There are many ways to implement elitism, but here we define elitism as the preservation of the best individual at each generation in a separate partition of the population space. This enlarges the population size by one individual; the elite individual increases the population size from N to N+1. However, note that the population size is still constant (i.e., equal to N+1) from one generation to the next. The elite individual does not take part in recombination or mutation but is maintained separately from the other N-members of the population. At each generation, if an individual in the N-member main population is better than the elite individual, then the elite individual is replaced with a copy of the better individual.Relative to a standardN-member BBO population, elite BBO increases the number of possible population distributions by a factor of n, which is the search space size. That is, each possible population distribution of the N-member main population could also include one of n elite individuals. The number of possible population distributions increases by a factor of n, from T (see Section 3.1) to nT. We order these new states so that each group of n states has the same elite individual. Also, the elite individual in the mth group of n states is the mth best individual in the search space, for m=1,…,n.The elitist-preserving process can be represented by an upgrade transition matrixO, which contains the probabilities that each population distribution of the (N+1)-member population transitions to some other population distribution after the elitist-preserving step. That is, the element in the ith row and jth column of O, denoted as O(i,j), is the probability that the ith population distribution transitions to the jth population distribution after the step in which the elite individual is replaced with the best individual from the N-member main population. The upgrade matrix is similar to the one in [29]; it does not include the effects of migration or mutation but only includes the elitism-preserving step. The upgrade matrix only includes the probability of changing the elite individual; it does not include the probability of changing the N-member main population, since it does not include migration or mutation. If there are no individuals in the N-member main population that are better than the elite individual, then the elite individual does not change. The structure of the upgrade matrix O can be written as
(22)O=(O110⋯0O21O22⋱⋮⋮⋱⋱0On1⋯On,n-1Onn),
where each Oij matrix is T×T, where T is the number of population distributions in an EA with a population size of N and search space cardinality of n. O11 is the identity matrix since the first n population distributions have the global optimum as their elite individual, and the elite individual can never be improved from the global optimum. Matrices Oaa with a≥2 are diagonal matrices composed of all zeros and ones. Since the population distributions are ordered by grouping common elite individuals and since elite individuals in the population distribution ordering are in order of decreasing fitness, the super block diagonals in O are zero matrices as shown in (22); that is, there is zero probability that the ith population distribution transitions to the jth population distribution if i<j. So the Markov chain of elite BBO can be described as
(23)P+=diag⁡(P)O=(P0⋯00P⋱⋮⋮⋱⋱00⋯0P)(O110⋯0O21O22⋱⋮⋮⋱⋱0On1⋯On,n-1Onn)=(PO110⋯0PO21PO22⋱⋮⋮⋱⋱0POn1⋯POn,n-1POnn),
where P is the T×T transition matrix described in Section 3.1.Example 12.
To explain the update matrixO described in (22), a simple example is presented. Suppose there exists a search space consisting of n=3 individuals which are x={x1,x2,x3} where the fitness of x1 is the lowest and the fitness of x3 is the highest. Suppose the main population size is N=1, so the elitist population size is N+1=2. So there are nine possible populations before the elitist-preserving step, and they are {C1,C2,C3,C4,C5,C6,C7,C8,C9} = {{x3,x3}, {x3,x2}, {x3,x1}, {x2,x3}, {x2,x2}, {x2,x1}, {x1,x3}, {x1,x2}, {x1,x1}}. Note that the first element in each population is the elite individual and the last N element (N=1 in this example) is the main population. Also note that the populations are ordered in such a way that the first three have the most fit individual as their elite individual, the next three have the second most fit individual as their elite individual, and the last three have the least fit individual as their elite individual. The update matrix O is a 9×9 matrix.The populationC1={x3,x3} transitions to the population C1={x3,x3} with probability 1; that is, O(1,1) = 1. Population C1 cannot transition to any other population Ci(i≠1); that is, O(1,i) = 0 for i≠1. Similarly, population C2={x3,x2} transitions to C2 with probability 1 since the elite x3 is better than the main-member population x2; therefore, O(2,2)=1, and O(2,i)=0 for i≠2. Continuing with this reasoning, we obtain the O matrix as follows:
(24)O=[O11O12O13O21O22O23O31O32O33]=[101001100000001000001100000000001000000000001],
where each Oij matrix is 3×3 and the blank elements are 0.Now we consider the convergence of thea*(t)-chain, which is the sequence of elite individuals in the elite BBO algorithm. If the elite individual is equal to the global optimum, we call this an absorbing state of the a*(t)-chain. Recall that the elite individual in elite BBO can only be replaced by one with better fitness. Therefore, the a*(t)-chain of elite BBO contains three classes of states: (1) at least one absorbing state, (2) nonabsorbing states which transition to absorbing states in one step, and (3) nonabsorbing states which transition to nonabsorbing states in one step. So the transition matrix P^ of the a*(t)-chain, which we introduced in (18)–(21), can be written as
(25)P^=(Ik0RQ),
where Ik is a k×k unit matrix corresponding to optimal individuals (k is the number of optima), R is a matrix of order (|S|-k)×k corresponding to nonabsorbing states that transition to absorbing states (|S| is the cardinality of the state space S, so |S|-k is the number of nonabsorbing states), and Q is a matrix of order (|S|-k)×(|S|-k) corresponding to nonabsorbing states that transition to nonabsorbing states. The matrix P^ of (25) has the same structure as the matrix P^ in Theorem 11. It follows from Theorem 11 that the a*(t)-chain of elite BBO is globally convergent.These results are similar to the canonical GA [29], which is proven to never converge to the global optimum, but elitist variants of which are proven to converge to the global optimum. We sum up these results in the following corollary.Corollary 13.
BBO with migration and mutation does not converge to any of the global optima, but elite BBO, which preserves the best individual at each generation, converges to the global optimum.Proof.
This is an immediate consequence of Theorems6 and 10 (the nonconvergence of BBO without elitism), Theorem 11 (the convergence of BBO with elitism), and the discussion above.


## 3.5. Convergence Rate

The previous subsection analyzed the convergence properties of elite BBO, and this subsection discusses its convergence rate. The transition matrix of elite BBO aftert steps can be found from (25) as follows:
(26)P^t=(Ik0NtRQt),
where Nt=I+Q+Q2+⋯+Qt-1=(I-Q)-1(I-Qt). If ∥Q∥<1 the limiting distribution of the Markov chain of BBO can be found from P^∞, which can be written as
(27)P^∞=lim⁡t→∞P^t=[Ik0(I-Q)-1R0].
Modified BBO with elitism has been proven to converge to a global optimum in the previous subsection, and there exists a limiting distribution π*=π(0)P^∞, where π(0) = [π1(0), …, πk(0), πk+1(0),…, π|s|(0)]=[π-1(0),π-2(0)] (recall that k is the number of global optima). The convergence rate estimate of elite BBO can be obtained as follows.Theorem 14.
If∥Q∥=ρ<1, the convergence rate of elite BBO satisfies ∥π(t)-π*∥≤O(ρt).Proof.
Consider(28)∥π(t)-π*∥=∥π(0)P^t-π(0)P^∞∥=∥π(0)[[Ik0NtRQt]-[Ik0(I-Q)-1R0]]∥=∥[π-1(0)π-2(0)][00(Nt-(I-Q)-1)RQt]∥=∥π-2(0)[(Nt-(I-Q)-1)RQt]∥≤∥π-2(0)∥(∥Nt-(I-Q)-1∥·∥R∥+∥Qt∥)=∥π-2(0)∥(∥(I-Q)-1(I-Qt)-(I-Q)-1∥mmmmmmm·∥R∥+∥Qt∥)=∥π-2(0)∥(∥(I-Q)-1-(I-Q)-1Qt-(I-Q)-1∥mmmmmmm·∥R∥+∥Qt∥)≤∥π-2(0)∥(∥-(I-Q)-1∥·∥R∥+1)∥Qt∥≤∥π-2(0)∥(∥R∥1-∥Q∥+1)∥Q∥t=O(ρt).
Note that elite BBO is guaranteed to converge to a global optimum regardless of the initial state. In addition, note that we can improve the convergence rate bound by decreasing the parameter ρ. That is, reducing the number of nonabsorbing states which transition to other nonabsorbing states can accelerate the convergence of elite BBO. In spite of differences between GAs and BBO [12], we see from Theorem 14 that the convergence rate of BBO with elitism is very similar to that of GAs [22, Theorem 5].


## 4. Simulation Results

Theorem14 gives the upper bound of the convergence rate estimate of elite BBO. In this section, we use simulation experiments to confirm this theorem. Note that in (28) the parameter ρ is a norm: ρ=∥Q∥. Here we define ∥·∥ as the infinity norm ∥·∥∞; that is, ∥Q∥∞=max⁡i⁡(∑j=1nqij), where qij is the element in the ith row and the jth column of matrix Q. Now note that the transition matrix P^ in (25) can be obtained from (11) and (23) using elementary matrix transformations. We can thus use Theorem 11 to check for BBO convergence, and we can use Theorem 14 to estimate the convergence rate of BBO. That is, we define π(t)-π* as the error between a BBO population distribution and a distribution that includes at least one optimal solution. We then define the convergence criteria as an arbitrarily small error (e.g., ∥π(t)-π*∥=10-6). We can then estimate the time t to convergence from (28) as follows:
(29)t≈logρ10-6.
Test functions are limited to three-bit problems with a search space cardinality of eight and a population size of four. The three fitness functions that we examine are
(30)f1=(12232334),f2=(42232334),f3=(41121223),
where f1 is a unimodal one-max problem, f2 is a multimodal problem, and f3 is a deceptive problem. Note that all three problems are to be maximized. Fitness values are listed in binary order, so the first element of each fitness function corresponds to the bit string 000, the second element corresponds to the bit string 001, and so on. For the BBO parameters, we use a maximum immigration rate and maximum emigration rate of 1, and we use linear migration curves as described in (2). We test elite BBO with three different mutation rates which are applied to each bit in each individual at each generation: 0.1, 0.01, and 0.001. Note that we do not test with a zero-mutation rate because the theory in this paper requires that the mutation rate be positive (see Theorem 5). Convergence is not guaranteed unless the mutation rate is positive.Numerical calculations show that the transition matrices for these three problems satisfy the convergence conditions of Theorem11, which indicates that the BBO algorithm converges to one or more of the global optima. As a heuristic test of Theorem 14, we use simulations to record the generation number of first obtaining a population in which all individuals are optimal, and all results are computed from 25 independent runs. Tables 1, 2, and 3 show comparisons of the theoretical convergence time t, the corresponding parameter ρ, and the generation number of first finding an all-optimal population, averaged over 25 independent simulations.Table 1
Convergence rate comparisons for the three-bit unimodal one-max problemf1.
Mutation rate
Theoretical analysis
Average generation number using simulation
Parameterρ
Convergence timet
BBO
GA
BBO
GA
BBO
GA
0.1
0.87
0.85
99.20
90.45
87.58
82.39
0.01
0.68
0.73
35.82
44.16
42.12
50.18
0.001
0.30
0.44
11.47
16.81
11.63
17.37
The table shows the convergence timet in seconds, the corresponding ρ based on (11) and (28), and the generation number of first finding an all-optimal population using BBO and GA, averaged over 25 Monte Carlo simulations.Table 2
Convergence rate comparison for the three-bit multimodal problemf2.
Mutation rate
Theoretical analysis
Average generation number using simulation
Parameterρ
Convergence timet
BBO
GA
BBO
GA
BBO
GA
0.1
0.61
0.59
27.95
25.16
26.09
24.87
0.01
0.11
0.17
6.26
7.95
6.17
7.72
0.001
0.04
0.05
4.29
4.76
4.18
4.63
The table shows the convergence timet in seconds, the corresponding ρ based on (11) and (28), and the generation number of first finding an all-optimal population using BBO and GA, averaged over 25 Monte Carlo simulations.Table 3
Convergence rate comparison for the three-bit deceptive problemf3.
Mutation rate
Theoretical analysis
Average generation number using simulation
Parameterρ
Convergence timet
BBO
GA
BBO
GA
BBO
GA
0.1
0.96
0.95
338.43
270.13
289.55
256.37
0.01
0.48
0.52
18.82
21.49
20.65
22.44
0.001
0.25
0.29
9.97
11.14
9.92
11.08
The table shows the convergence timet in seconds, the corresponding ρ based on (11) and (28), and the generation number of first finding an all-optimal population using BBO and GA, averaged over 25 Monte Carlo simulations.Tables1–3 show time to convergence and time to finding an optimum, for both BBO and GA. The table confirms the statement following Theorem 14 that the convergence behavior of BBO is similar to that of GA. The tables show that GA converges slightly faster than BBO for high mutation rates, but BBO converges slightly faster for low mutations, and this latter behavior is more important in practice because low mutations rates provide faster convergence.Several things are notable about the results in Tables1–3. First, the mutation rate affects the convergence rate of BBO. For all test problems, the convergence rate improves when the mutation rate decreases. We can accelerate the convergence of BBO by decreasing the mutation rate. This may provide practical guidance for BBO tuning for real-world problems. Second, by analyzing the relationship of the parameter ρ and the convergence time t in Tables 1–3, we see that the convergence time t is exponentially related to the parameter ρ, as predicted by Theorem 14. Third, the theoretical results and simulation results match well for most of the test problems, which confirms the convergence rate estimate provided by Theorem 14.The three-bit problems analyzed above are small, and the results do not tell us how convergence rates scale with the problem dimension. Also, the transition matrix grows faster than exponentially with the problem dimension and the population size [11], so realistically sized problems cannot be directly analyzed with the methods in this paper. However, our methods could be used to study the effect of BBO tuning parameters on small problems, which could provide guidance for larger, real-world problems. Also, similar population distributions could be grouped into the same Markov model state to reduce the transition matrix dimension of large problems to a manageable size [16], which could make the methods in this paper practical for realistically sized problems.


## 5. Conclusion

In this paper we modeled BBO as a homogeneous finite Markov chain to study convergence, and we obtained new theoretical results for BBO. The analysis revealed that BBO with only migration and mutation does not converge to the global optimum. However, an elite version of BBO that maintains the best solution in the population from one generation to the next converges to a population subset containing at least one globally optimal solution. In other words, BBO with elitism will converge to the global optimum in any binary optimization problem.In addition, an upper bound for the BBO convergence rate was obtained in Theorem14. We used simulations to confirm this theorem for a unimodal one-max problem, a multimodal problem, and a deceptive problem. The results in this paper are similar to those of the canonical GA [29], and so our results are not surprising, but this paper represents the first time that such results have been formalized for BBO.The results in this paper are limited to binary problems but can easily be extended to discrete problems with any finite alphabet. This is due to the simple fact that any discrete problem can be reduced to a binary problem.For future work there are several important directions. First, it is of interest to study how to improve BBO convergence time and robustness based on these results. The second important direction for future work is to study the asymptotic convergence of variations of BBO, including partial emigration-based BBO, total immigration-based BBO, and total emigration-based BBO [12]. The theorems in this paper provide the foundation to study these variations, so we do not need additional theoretical tools to analyze their convergence. The third important direction for future work is to develop hybrid BBO, which combines BBO with other EAs, and study their convergence behaviors using the theory presented here. Finally, it would be of interest to extend these results to continuous optimization problems, which are the types of problems to which real-world implementations of BBO are typically applied.



---

*Source: 147457-2014-05-22.xml*