演算法背景
人工蜂群演算法 (Artificial Bee Colony, ABC) 是由 Karaboga 於 2005 年提出的一種新穎的基於叢集智慧的全域性優化演算法,其直觀背景來源於蜂群的採蜜行為。它的主要特點是不需要了解問題的特殊資訊,只需要對問題進行優劣的比較,通過各人工蜂個體的區域性尋優行為,最終在群體中使全域性最優值突現出來,有著較快的收斂速度。
蜜蜂是一種群居昆蟲,雖然單個昆蟲的行為極其簡單,但是由單個簡單的個體所組成的群體卻表現出極其複雜的行為。真實的蜜蜂種群能夠在任何環境下,以極高的效率從食物源(花朵)中採集花蜜;同時,它們能適應環境的改變。
搜尋流程
演算法的呼叫過程如下:
初始化所有蜜源
記錄最優蜜源
while:
僱傭蜂對所有蜜源進行鄰域搜尋(避免飢餓效應)
計算輪盤度,判斷蜜源質量
觀察蜂對優質蜜源進行鄰域搜尋(加速演算法收斂)
記錄最優蜜源
偵查蜂放棄枯竭蜜源進行全域性搜尋(跳出區域性最優)
記錄最優蜜源
end
其中僱傭蜂和觀察蜂有著相似的邏輯,特別在對指定蜜源進行鄰域搜尋時,兩者的邏輯是完全的一樣的:
- 基於原有蜜源進行鄰域突變
- 保證鄰域突變的有效性
- 若為優質突變,則進行蜜源替換
- 若為劣質突變,則進行蜜源開採
但是演算法的設計者們卻特意區分出兩種不同的邏輯,其原因可以從實現程式碼中看出。在進行領域搜尋時,對指定蜜源的選擇和限定是關鍵所在,它暗示了僱傭蜂和觀察蜂的區別以及承擔的不同角色。
首先對於僱傭蜂的角色,其指定蜜源的方式簡單粗暴,對每一個蜜源進行遍歷指定。通過這種方式進行鄰域搜尋,是建立整個演算法的基礎核心。
而對於觀察蜂角色,它是根據輪盤賭策略進行蜜源的指定,也就是說,蜜源越是優質,其被指定的、被進行領域搜尋的概率就越高。通過這種正向反饋的過程,加速了整個演算法的收斂性,可以幫助我們在多個區域性中快速找到最優解。
如此看來觀察蜂似乎是僱傭蜂的進化版,觀察蜂似乎可以完全替代僱傭蜂?其實不然。觀察蜂角色在進行快速收斂、對優質蜜源進行了較多照顧的同時,劣質的蜜源可能會被忽略,從而產生飢餓效應。而僱傭蜂的角色讓所有的蜜源,不論優劣,都有能夠持續搜尋的機會,僱傭蜂和觀察蜂相互配合,使得演算法能夠均衡高效的執行。
最後的角色是偵查蜂。當持續的搜尋使得蜜源枯竭的時候,偵查蜂對枯竭的蜜源進行放棄,跳出原有的區域性空間,在全域性空間中隨機探索新的蜜源,之後轉變為僱傭蜂重複上述過程。
Java 實現
主函式
public class test {
// 建立蜂群
static beeColony bee = new beeColony();
public static void main(String[] args) {
int iter = 0;
int run = 0;
int j = 0;
double mean = 0;
//srand(time(NULL));
// 重複計算,獲得平均情況下的全域性最優解
for (run = 0; run < bee.runtime; run++) {
// 初始化蜜源
bee.initial();
bee.MemorizeBestSource();
for (iter = 0; iter < bee.maxCycle; iter++) {
// 僱傭蜂的邏輯: 更新所有蜜源(確保能跳出區域性蜜源)
bee.SendEmployedBees();
bee.CalculateProbabilities();
// 觀察蜂的邏輯: 好蜜源正向反饋(加速收斂的作用)
bee.SendOnlookerBees();
bee.MemorizeBestSource();
// 僅僅重製枯竭的蜜源 (區域性優質蜜源枯竭的快)
bee.SendScoutBees();
}
// 列印最優的座標,以及最終的函式值
for (j = 0; j < bee.D; j++) {
//System.out.println("GlobalParam[%d]: %f\n",j+1,GlobalParams[j]);
System.out.println("GlobalParam[" + (j + 1) + "]:" + bee.GlobalParams[j]);
}
//System.out.println("%d. run: %e \n",run+1,GlobalMin);
System.out.println((run + 1) + ".run:" + bee.GlobalMin);
bee.GlobalMins[run] = bee.GlobalMin;
mean = mean + bee.GlobalMin;
}
mean = mean / bee.runtime;
//System.out.println("Means of %d runs: %e\n",runtime,mean);
System.out.println("Means of " + bee.runtime + "runs: " + mean);
}
}
演算法函式
import java.lang.Math;
public class beeColony {
/* Control Parameters of ABC algorithm*/
int NP = 20; /* The number of colony size (employed bees+onlooker bees)*/
int FoodNumber = NP / 2; /*The number of food sources equals the half of the colony size*/
int limit = 100; /*A food source which could not be improved through "limit" trials is abandoned by its employed bee*/
int maxCycle = 2500; /*The number of cycles for foraging {a stopping criteria}*/
/* Problem specific variables*/
int D = 100; /*The number of parameters of the problem to be optimized*/
double lb = -5.12; /*lower bound of the parameters. */
double ub = 5.12; /*upper bound of the parameters. lb and ub can be defined as arrays for the problems of which parameters have different bounds*/
int runtime = 30; /*Algorithm can be run many times in order to see its robustness*/
int dizi1[] = new int[10];
double Foods[][] = new double[FoodNumber][D]; /*Foods is the population of food sources. Each row of Foods matrix is a vector holding D parameters to be optimized. The number of rows of Foods matrix equals to the FoodNumber*/
double f[] = new double[FoodNumber]; /*f is a vector holding objective function values associated with food sources */
double fitness[] = new double[FoodNumber]; /*fitness is a vector holding fitness (quality) values associated with food sources*/
double trial[] = new double[FoodNumber]; /*trial is a vector holding trial numbers through which solutions can not be improved*/
double prob[] = new double[FoodNumber]; /*prob is a vector holding probabilities of food sources (solutions) to be chosen*/
double solution[] = new double[D]; /*New solution (neighbour) produced by v_{ij}=x_{ij}+\phi_{ij}*(x_{kj}-x_{ij}) j is a randomly chosen parameter and k is a randomlu chosen solution different from i*/
double ObjValSol; /*Objective function value of new solution*/
double FitnessSol; /*Fitness value of new solution*/
int neighbour, param2change; /*param2change corrresponds to j, neighbour corresponds to k in equation v_{ij}=x_{ij}+\phi_{ij}*(x_{kj}-x_{ij})*/
double GlobalMin; /*Optimum solution obtained by ABC algorithm*/
double GlobalParams[] = new double[D]; /*Parameters of the optimum solution*/
double GlobalMins[] = new double[runtime];
/*GlobalMins holds the GlobalMin of each run in multiple runs*/
double r; /*a random number in the range [0,1)*/
/*a function pointer returning double and taking a D-dimensional array as argument */
/*If your function takes additional arguments then change function pointer definition and lines calling "...=function(solution);" in the code*/
// typedef double (*FunctionCallback)(double sol[D]);
/*benchmark functions */
// double sphere(double sol[D]);
// double Rosenbrock(double sol[D]);
// double Griewank(double sol[D]);
// double Rastrigin(double sol[D]);
/*Write your own objective function name instead of sphere*/
// FunctionCallback function = &sphere;
/*Fitness function*/
double CalculateFitness(double fun) {
double result = 0;
if (fun >= 0) {
result = 1 / (fun + 1);
} else {
result = 1 + Math.abs(fun);
}
return result;
}
/*The best food source is memorized*/
void MemorizeBestSource() {
int i, j;
for (i = 0; i < FoodNumber; i++) {
if (f[i] < GlobalMin) {
GlobalMin = f[i];
for (j = 0; j < D; j++)
GlobalParams[j] = Foods[i][j];
}
}
}
/*Variables are initialized in the range [lb,ub]. If each parameter has different range, use arrays lb[j], ub[j] instead of lb and ub */
/* Counters of food sources are also initialized in this function*/
void init(int index) {
int j;
for (j = 0; j < D; j++) {
r = ((double) Math.random() * 32767 / ((double) 32767 + (double) (1)));
Foods[index][j] = r * (ub - lb) + lb;
solution[j] = Foods[index][j];
}
f[index] = calculateFunction(solution);
fitness[index] = CalculateFitness(f[index]);
trial[index] = 0;
}
/*All food sources are initialized */
void initial() {
int i;
for (i = 0; i < FoodNumber; i++) {
init(i);
}
GlobalMin = f[0];
for (i = 0; i < D; i++)
GlobalParams[i] = Foods[0][i];
}
void SendEmployedBees() {
int i, j;
/*Employed Bee Phase*/
for (i = 0; i < FoodNumber; i++) {
/*The parameter to be changed is determined randomly*/
r = ((double) Math.random() * 32767 / ((double) (32767) + (double) (1)));
param2change = (int) (r * D);
/*A randomly chosen solution is used in producing a mutant solution of the solution i*/
r = ((double) Math.random() * 32767 / ((double) (32767) + (double) (1)));
neighbour = (int) (r * FoodNumber);
/*Randomly selected solution must be different from the solution i*/
// while(neighbour==i)
// {
// r = ( (double)Math.random()*32767 / ((double)(32767)+(double)(1)) );
// neighbour=(int)(r*FoodNumber);
// }
for (j = 0; j < D; j++)
solution[j] = Foods[i][j];
/*v_{ij}=x_{ij}+\phi_{ij}*(x_{kj}-x_{ij}) */
r = ((double) Math.random() * 32767 / ((double) (32767) + (double) (1)));
solution[param2change] = Foods[i][param2change] + (Foods[i][param2change] - Foods[neighbour][param2change]) * (r - 0.5) * 2;
/*if generated parameter value is out of boundaries, it is shifted onto the boundaries*/
if (solution[param2change] < lb)
solution[param2change] = lb;
if (solution[param2change] > ub)
solution[param2change] = ub;
ObjValSol = calculateFunction(solution);
FitnessSol = CalculateFitness(ObjValSol);
/*a greedy selection is applied between the current solution i and its mutant*/
if (FitnessSol > fitness[i]) {
/*If the mutant solution is better than the current solution i, replace the solution with the mutant and reset the trial counter of solution i*/
trial[i] = 0;
for (j = 0; j < D; j++)
Foods[i][j] = solution[j];
f[i] = ObjValSol;
fitness[i] = FitnessSol;
} else { /*if the solution i can not be improved, increase its trial counter*/
trial[i] = trial[i] + 1;
}
}
/*end of employed bee phase*/
}
/* A food source is chosen with the probability which is proportioal to its quality*/
/*Different schemes can be used to calculate the probability values*/
/*For example prob(i)=fitness(i)/sum(fitness)*/
/*or in a way used in the metot below prob(i)=a*fitness(i)/max(fitness)+b*/
/*probability values are calculated by using fitness values and normalized by dividing maximum fitness value*/
void CalculateProbabilities() {
int i;
double maxfit;
maxfit = fitness[0];
for (i = 1; i < FoodNumber; i++) {
if (fitness[i] > maxfit)
maxfit = fitness[i];
}
for (i = 0; i < FoodNumber; i++) {
prob[i] = (0.9 * (fitness[i] / maxfit)) + 0.1;
}
}
// just add choose tactic, the center logic is as same as sendEmployedBees
void SendOnlookerBees() {
int i, j, t;
i = 0;
t = 0;
/*onlooker Bee Phase*/
while (t < FoodNumber) {
r = ((double) Math.random() * 32767 / ((double) (32767) + (double) (1)));
if (r < prob[i]) /*choose a food source depending on its probability to be chosen*/ {
t++;
/*The parameter to be changed is determined randomly*/
r = ((double) Math.random() * 32767 / ((double) (32767) + (double) (1)));
param2change = (int) (r * D);
/*A randomly chosen solution is used in producing a mutant solution of the solution i*/
r = ((double) Math.random() * 32767 / ((double) (32767) + (double) (1)));
neighbour = (int) (r * FoodNumber);
/*Randomly selected solution must be different from the solution i*/
while (neighbour == i) {
//System.out.println(Math.random()*32767+" "+32767);
r = ((double) Math.random() * 32767 / ((double) (32767) + (double) (1)));
neighbour = (int) (r * FoodNumber);
}
for (j = 0; j < D; j++)
solution[j] = Foods[i][j];
/*v_{ij}=x_{ij}+\phi_{ij}*(x_{kj}-x_{ij}) */
r = ((double) Math.random() * 32767 / ((double) (32767) + (double) (1)));
solution[param2change] = Foods[i][param2change] + (Foods[i][param2change] - Foods[neighbour][param2change]) * (r - 0.5) * 2;
/*if generated parameter value is out of boundaries, it is shifted onto the boundaries*/
if (solution[param2change] < lb)
solution[param2change] = lb;
if (solution[param2change] > ub)
solution[param2change] = ub;
ObjValSol = calculateFunction(solution);
FitnessSol = CalculateFitness(ObjValSol);
/*a greedy selection is applied between the current solution i and its mutant*/
if (FitnessSol > fitness[i]) {
/*If the mutant solution is better than the current solution i, replace the solution with the mutant and reset the trial counter of solution i*/
trial[i] = 0;
for (j = 0; j < D; j++)
Foods[i][j] = solution[j];
f[i] = ObjValSol;
fitness[i] = FitnessSol;
} else { /*if the solution i can not be improved, increase its trial counter*/
trial[i] = trial[i] + 1;
}
} /*if */
i++;
if (i == FoodNumber)
i = 0;
}/*while*/
/*end of onlooker bee phase */
}
/*determine the food sources whose trial counter exceeds the "limit" value. In Basic ABC, only one scout is allowed to occur in each cycle*/
void SendScoutBees() {
int maxtrialindex, i;
maxtrialindex = 0;
// 遍歷當前所有食物源,尋找最大的trail
for (i = 1; i < FoodNumber; i++) {
if (trial[i] > trial[maxtrialindex])
maxtrialindex = i;
}
if (trial[maxtrialindex] >= limit) {
init(maxtrialindex);
}
}
// target function for calculate
double calculateFunction(double sol[]) {
return Rastrigin(sol);
}
double sphere(double sol[]) {
int j;
double top = 0;
for (j = 0; j < D; j++) {
top = top + sol[j] * sol[j];
}
return top;
}
double Rosenbrock(double sol[]) {
int j;
double top = 0;
for (j = 0; j < D - 1; j++) {
top = top + 100 * Math.pow((sol[j + 1] - Math.pow((sol[j]), (double) 2)), (double) 2) + Math.pow((sol[j] - 1), (double) 2);
}
return top;
}
double Griewank(double sol[]) {
int j;
double top1, top2, top;
top = 0;
top1 = 0;
top2 = 1;
for (j = 0; j < D; j++) {
top1 = top1 + Math.pow((sol[j]), (double) 2);
top2 = top2 * Math.cos((((sol[j]) / Math.sqrt((double) (j + 1))) * Math.PI) / 180);
}
top = (1 / (double) 4000) * top1 - top2 + 1;
return top;
}
double Rastrigin(double sol[]) {
int j;
double top = 0;
for (j = 0; j < D; j++) {
top = top + (Math.pow(sol[j], (double) 2) - 10 * Math.cos(2 * Math.PI * sol[j]) + 10);
}
return top;
}
}
參考連結
[1] https://abc.erciyes.edu.tr/
[2] https://en.wikipedia.org/wiki/Test_functions_for_optimization
[3] https://www.pianshen.com/article/729179041/
[4] https://www.jianshu.com/p/ebd436d27cf8