электронное научно-техническое издание
НАУКА и ОБРАЗОВАНИЕ
Эя № ФС 77 - 30569. Государственная регистрация №0420900025. ISSN 1994-0408
Исследование эффективности метода пчелиного роя в задаче глобальной оптимизации
# 08, август 2010
авторы: Гришин А. А., Карпенко А. П.
УДК 519.6
МГТУ им. Н.Э. Баумана, grishinalexandera@gmail. com;
Введение
Среди задач непрерывной конечномерной оптимизации самым важным с практической точки зрения и, одновременно, самым сложным является класс задач глобальной оптимизации. Методы решения задачи глобальной оптимизации делятся на детерминированные стохастические и эвристические методы [1].
Эвристические методы являются относительно новым и быстро развивающимся классом методов глобальной оптимизации. Среди этих методов выделяются эволюционные и поведенческие (имитационные) методы.
Поведенческие методы относятся к мультиагентным методам, основанным на моделировании интеллектуального поведения колоний агентов (Swarm Intelligence) [2, 3]. В природе таким интеллектом обладают группы общественных насекомых, например, колонии термитов, муравьёв, пчёл, некоторых видов ос.
Динамика популяции общественных насекомых определяется взаимодействиями насекомых друг с другом, а также с окружающей средой. Эти взаимодействия осуществляются посредством различных химических и/или физических сигналов, например, феромонов, выделяемых муравьями.
Метод пчелиного роя (Bees Algorithm) является одним из новейших методов, относящихся к рассматриваемому направлению. Первые статьи, в которых был предложен данный метод, были опубликованы в 2005 г. [4, 5]. Метод представляет собой эвристический итеративный мультиагентный метод случайного поиска, основная идея которого состоит в моделировании поведения пчёл при поиске нектара.
Среди недостатков метода пчелиного роя упоминания заслуживает большое число свободных параметров метода, от значений которых, с одной стороны, зачастую сильно зависит эффективность метода, а, с другой стороны, отсутствуют какие-либо содержательные основания для выбора этих значений.
Работа посвящена исследованию эффективности варианта метода пчелиного роя, предложенного в публикации [6]. Приведены бионические предпосылки метода, представлена схема используемого варианта метода. Программная реализация метода выполнена в среде программирования Python, в котором для ускорения вычислений над векторами и матрицами использованы модули расширения Psyco и Numpy [7]. Тестирование разработанного программного обеспечения выполнено на двумерной функции Шекеля. Исследование эффективности метода и программного обеспечения выполнено на трех тестовых функциях из известного пакета CES (Congress of Evolutionary Computing) [8]. Указанные функции являются представителями трех классов функций: одноэкстремальные овражные функции; многоэкстремальные функции, имеющие небольшое число экстремумов; многоэкстремальные функции, имеющие большое число экстремумов.
Отметим, что известны примеры успешного применения метода пчелиного роя для решения ряда прикладных задач: задачи календарного планирования [9]; задачи коммивояжёра [10]; транспортной задачи [11] и др. [12, 13].
1. Бионические предпосылки
Для описания поведения пчёл в природе используются три следующих основных понятия [14].
Источник нектара характеризуется своей полезностью, которая определяется такими факторами, как удалённость от улья, концентрация нектара, удобство его добычи.
Занятые фуражиры - пчелы, которые «связаны» с одним из источников нектара, т.е. добывают на нем нектар. Занятые фуражиры владеют следующей информацией о «своем» источнике нектара: направление от улья на источник; полезность источника.
Незанятые фуражиры - пчелы-разведчики, которые осуществляют поиск источников нектара для их использования, а также пчелы-наблюдатели, которые в данное время выполняют некоторые работы в улье.
Каждый незанятый фуражир может полететь к источнику нектара, следуя за пчелой-разведчиком, которая нашла путь к такому источнику. Пчела-разведчик выполняет «вербовку» незанятых пчел с помощью танца на специальной площадке улья - области танцев. Завербованная пчела следует за соответствующей пчелой-разведчиком к области с нектаром и становится, таким образом, занятым фуражиром.
Занятый фуражир после добычи нектара возвращается в улей и оставляет его там. После этого данный фуражир может выполнить одно из следующих действий: оставить «свой» источник нектара и стать незанятым фуражиром; продолжить заготовку нектара из прежнего источника, не вербуя других пчел; выполнить вербовку. Пчела выбирает одно из указанных действий по некоторому вероятностному закону.
Одновременно в пределах области танцев разные пчелы могут "рекламировать" различные источники нектара. Механизмы принятия решений, в соответствии с которыми пчела решает следовать за той или иной пчелой-вербовщиком, исследованы недостаточно хорошо. Логично предположить, что вероятность вербовки тем или иным образом определяется полезностью соответствующего источника нектара [14].
Таким образом, самоорганизация пчелиного роя основывается на четырёх следующих основных механизмах.
1). Положительная обратная связь - на основе информации, полученной от других пчел, пчела летит к одному из источников нектара.
2). Отрицательная обратная связь - основываясь на информации, полученной от других пчёл, данная пчела может решить, что «ее» источник нектара значительно хуже других найденных источников, и оставить этот источник.
3). Случайность - вероятностный поиск пчёлами-разведчиками новых источников нектара.
4). Множественность взаимодействия - информация об источнике нектара, найденном одной пчелой, передается многим другим пчелам улья.
2. Постановка задачи и схема используемого метода
Рассмотрим задачу глобальной условной оптимизации
maxF(X) = F(X*) = F*, X е D с RK,
X
где X = (x1,x2,...,xK) - вектор варьируемых параметров,
D = {X | xj < xt < x+,i е [1: k]} -множество допустимых значений этого вектора.
Обозначим пчелиный рой B = [Ba,a е [1: Z]}, где Ba - пчела (агент), Z -число агентов в рое. Положение пчелы Ba в момент времени t = 0,1,...
определяется вектором ее координат Ха г = (ха г 1,ха г1,..,хагк). Пусть
Бр е В, Р е [1: £ ] - пчела-разведчик; £ < 2.
Схема используемого варианта метода роя пчел имеет следующий вид. На первом шаге метода в точки со случайными координатами X р 0 е D,
отправляются пчелы-разведчики. В зависимости от значений целевой функции F(X) в этих точках, в области D выделяются два вида участков
(подобластей) dр:
• п лучших участков, которые соответствуют наибольшим значениям
целевой функции;
• т перспективных участков, соответствующих значениям целевой
функции, наиболее близким к наилучшим значениям. Подобласть dр, Р е [1: £] называется подобластью локального поиска
и представляет собой гиперкуб в пространстве Rк с центром в точке Xр 0
и длинами сторон, равными 2А. Здесь А - параметр, называемый размером области локального поиска.
Если оказалось, что евклидово расстояние
Хр,0 Ху,0
между двумя
агентами-разведчиками Бр,Б* е В, р ф Л меньше некоторой фиксированной
величины, то, вообще говоря, возможны два следующих варианта метода, из которых в работе используется второй вариант:
• поставить в соответствие этим агентам два различных
пересекающихся участка dр, dr (лучших и/или перспективных);
• поставить в соответствие тем же агентам один участок, центр
которого находится в точке, соответствующей агенту с большим значением целевой функции. В каждый из лучших и перспективных участков посылается по N и по М агентов, соответственно. Координаты этих агентов в указанных участках определяются случайным образом.
Отметим, что вариантом рассмотренного решения является посылка в указанные подобласти не фиксированных чисел агентов, а чисел, пропорциональных соответствующим значениям целевой функции. Размеры подобластей, в которые посылаются агенты, могут уменьшаться с ростом числа итераций с тем, чтобы в каждой из подобластей решение сходилось к локальному максимуму целевой функции в этой подобласти.
На основе анализа значений целевой функции, соответствующих всем агентам роя, после некоторого числа итераций находятся n новых лучших и m новых перспективных участков. В качестве критерия окончания итераций можно использовать достижение заданного количества итераций T. Можно также заканчивать итерации, если в течение т < T итераций не удалось увеличить максимальное достигнутое значение целевой функции. Здесь т - параметр останова.
Более сложный вариант метода пчелиного роя рассмотрен, например, в работе [15].
3. Программная реализация метода
Программная реализация метода пчелиного роя выполнена на языке программирования Python, который представляет собой высокоуровневый интерпретируемый язык программирования общего назначения. Синтаксис ядра Python минималистичен, в то же время, его стандартная библиотека включает в себя большое число функций [7].
Python поддерживает несколько парадигм программирования, в том числе структурное, объектно-ориентированное, функциональное, императивное и аспектно-ориентированное. Основные архитектурные черты языка — динамическая типизация, автоматическое управление памятью, полная интроспекция, механизм обработки исключений, поддержка многопоточных вычислений и удобные высокоуровневые структуры данных. Программный код на языке Python организуется в функции и классы,
которые могут объединяться в модули, которые в свою очередь могут быть объединены в пакеты.
Эффективность Python позволяют значительно повысить его модули расширения Psyco и Numpy [7]. Во время исполнения программы на языке Python, модуль расширения Psyco может выборочно заменять части интерпретируемого байткода Python оптимизированным машинным кодом.
Поскольку Python - интерпретируемый язык, математические алгоритмы часто работают в нём гораздо медленнее, чем в компилируемых языках программирования, таких как C или даже Java. Модуль расширения Numpy решает эту проблему для большого числа вычислительных алгоритмов, обеспечивая поддержку многомерных массивов, а также функций и операторов для работы с ними. В результате, любой алгоритм который может быть выражен, в основном, как последовательность операций над векторами и матрицами, выполняет также быстро, как эквивалентный код написанный на языке C.
Программа Bee, реализующая метод роя пчел, содержит следующие основные модули: pybee.py; beetest.py; beeexamples.py; beetestfunc.py.
Модуль pybee включает в себя следующие классы: floatbee - базовый класс; hive - класс, реализующий основные операции метода; statistic -вспомогательный класс, предназначенный для накопления результатов итераций метода.
Модуль beetest предназначен для задания свободных параметров метода.
С помощью модуля beeexamples задаются границы множества допустимых значений вектора варьируемых параметров D, а также границы областей локального поиска.
Модуль beetestfunc реализует вспомогательные функции, обеспечивающие визуализацию результатов вычислений.
4. Тестирование разработанного программного обеспечения
Для тестирования метода и разработанного программного обеспечения использовалась двумерная функция Шекеля (Shekel) [8]
F (х, у) =-^-2 +-1-2 +-1-2
1 + (х - 2)2 + (y -10)2 2 + (х -10)2 + (y -15)2 2 + (х -18)2 + (у - 4)2
Поиск выполнялся в области D = {(х,у)| х е[0,20],у е[0,20]}. Здесь и далее для компактности записи варьируемые параметры х1, х2 заменены параметрами х, у . Вид функции иллюстрирует рисунок 1.
Рис. 1 - Поверхность двумерной функции Шекеля F (х, у)
Легко видеть, что глобальный максимум в рассматриваемой задаче достигается в точке с координатами (2,10) и значение функции в этой точке
равно F = 1,01439037.
Использовались следующие значения свободных параметров метода:
• размер области локальный поиск, Д = 1;
• число лучших участков п = 2;
• число перспективных участков т = 3;
• число агентов, отправляемых на лучшие участки N = 7; Электронный журнал, №8 Август 2010г. http://technomag.edu.ru/ Страница 8
• число агентов, отправляемых на перспективные участки М = 4;
• максимальное количество итераций Т = 100.
Расположение агентов после выполнения первой, четвертой и 17-ой итераций иллюстрируют рисунки 2 - 4.
У -1-\-1-1-
з-о I
1Л
Рис. 2 - Распределение агентов: t = 1
у
о
ю ю х
Рис. 3 - Распределение агентов: t = 4
Рисунки показывают, что, как и следовало ожидать, с ростом номера итераций агенты все в большей степени концентрируются вблизи глобального максимума целевой функции.
На рисунке 5 представлена зависимость найденного программой максимального значения целевой функции от номера итерации t. Рисунок
показывает, что уже поле примерно 20 итераций программа с высокой точностью находит глобальный максимум.
Рис. 4 - Распределение агентов: t = 17
Рис. 5 - Максимальное достигнутое значение целевой функции F (x, y)
Таким образом, тестирование программы Bee, реализующей метод пчелиного роя, подтвердило работоспособность используемых метода, алгоритма и программного обеспечения.
x
5. Исследование эффективности метода
Исследование эффективности метода и программного обеспечения выполнено на трех тестовых функциях из пакета тестовых функций CES (Congress of Evolutionary Computing) - функции Розенброка, функции Химмельблау и функции Растригина [8].
Если не оговорено противное, имеются в виду следующие значения свободных параметров метода: число агентов-разведчиков S = 16; размер области локального поиска А = 0,2 ; число лучших участков n = 2 ; число перспективных участков m = 3 ; число агентов, отправляемые на лучшие участки N = 7 ; число агентов, отправляемых на перспективные участки M = 4; максимальное число итераций T = 500 ; параметр останова т = 20.
В процессе исследования варьировались значения следующих параметров:
• число агентов-разведчиков S ;
• размер области локального поиска А ;
• параметр останова т .
В качестве критерия окончания итераций во всех случаях использовано неулучшение решения в течение т последовательных итераций. Число итераций, при котором в соответствие с этим критерием, завершается итерационный процесс, обозначается te.
Эффективность метода в значительной мере зависит от начального расположения агентов-разведчиков. Чтобы избавиться от этой зависимости, в каждом исследовании производилось по 30 запусков программы, отличающихся начальным положением агентов. Приведенные ниже результаты получены путем усреднения исследуемых характеристик метода по этим запускам.
5.1. Функция Розенброка. Рассмотрена двумерная функция, обратная функции Розенброка (рисунок 6)
F(x,y) = -(100(y - x2)2 + (1 - x)2).
X . П
X
а) б)
Рис.6 -Поверхность (а) и линии уровня (б) двумерной функции Розенброка
Глобальный максимум этой функции достигается в точке с координатами (1,1) и равен нулю; область поиска представляет собой квадрат В = {(х,у) | х е [-5,5],у е [-5,5]}.
Варьирование числа агентов-разведчиков. Число агентов-разведчиков £ последовательно принималось равным 1, 2, 4, 8, 16, 32. Результаты исследования иллюстрируют рисунки 7, 8.
1 2 4 8 16 32
Число агентов-разведчиков
Рис. 7 - Функция Розенброка: число итераций te в функции числа агентов-
разведчиков £
1 2 4 8 16 32
Число агентов-разведчиков
Рис. 8 - Функция Розенброка: погрешность решения е, как функции числа
агентов-разведчиков 5
Рисунки показывают, что для функции Розенброка при увеличении числа агентов-разведчиков число итераций te и точность найденного решения е уменьшаются. Здесь и далее под точностью решения понимается его абсолютной погрешностью.
Варьирование размера области локального поиска. Исследование выполнено при значениях параметра Д, равных 0,1; 0,2; 0,5; 1,0. Результаты исследования иллюстрируют рисунки 9, 10.
0,1 0,2 0,5 1,0
Размер области локального поиска
Рис. 9 - Функция Розенброка: число итераций te в функции размера области
локального поиска Д
Рис. 10 - Функция Розенброка: погрешность решения s, как функция размера области локального поиска А
Рисунки 9, 10 показывают, что для функции Розенброка при уменьшении размеров области локального поиска, точность найденного решения s повышается, однако при этом увеличивается число итераций te.
Варьирование параметра останова. Число итераций т, за которое найденное решение не удалось улучшить, последовательно принималось равным 5, 10, 20, 50. Полученные результаты иллюстрируют рисунки 11, 12.
Рис. 11 - Функция Розенброка: число итераций te в функции параметра
останова т
Рис. 12 - Функция Розенброка: погрешность решения s, как функция
параметра останова г
Рисунки 11, 12 показывают, что для функции Розенброка с ростом величины г погрешность решения s уменьшается, а число итераций te возрастает.
5.2. Функция Химмельблау (Himmelblau). Рассмотрена также двумерная функция, обратная функции Химмельблау (рисунок 13)
F(х,у) = -((х2 + y -11)2 + (х + y2 - 7)2 ).
а)
б)
Рис. 13 - Поверхность (а) и линии уровня (б) двумерной функции
Химмельблау
Функция имеет четыре локальных максимума в точках (-2,805118; 3,131312), (-3,779310; -3,283186), (3,584428; -1,848126), (3,0; 2,0) и принимает
во всех этих точках нулевые значения. Область поиска принята равной D = {(х,у)| х е [-10,10],у е [-10,10]}. Если не оговорено противное, размер области локального поиска полагается равным Д = 0,4.
Варьирование числа агентов-разведчиков. Как и в п. 5.1, число агентов-разведчиков 5 принималось равным 1, 2, 4, 8, 16, 32. Соответствующие результаты представлены на рисунках 14, 15.
Рис. 14 - Функция Химмельблау: число итераций te в функции числа
агентов-разведчиков 5
Рисунок 15 - Функция Химмельблау: погрешность решения е, как функция
числа агентов-разведчиков 5
Рисунки 14, 15 показывают, что для функции Химмелблау число итераций te практически не зависит от числа агентов-разведчиков;
погрешность найденного решения е уменьшается с ростом этого числа.
Варьирование размера области локального поиска. Исследование выполнено для параметра А, равного 0,2; 0,4; 1,0; 2,0. Соответствующие результаты представлены на рисунках 16, 17.
Рис. 16 - Функция Химмельблау: число итераций te в функции размера
области локального поиска А
Рис. 17 - Функция Химмельблау: погрешности решения е, как функция размера области локального поиска А
Рисунки 16, 17 показывают, что для функции Химмельблау при увеличении размера области локального поиска Д погрешность полученного решения е растет, а число итераций te уменьшается.
Варьирование параметра останова т. Число итераций т, в течение которых решение не удалось улучшить, принимало в исследовании значения, равные 5, 10, 20, 50. Результаты исследования иллюстрируют рисунки 18, 19.
Рис. 18 - Функция Химмельблау: число итераций te в функции
параметра останова т
0,ОйО 0,04 5 0,040 0,035 О.ОЗО 0,025 0,020 0,015 О.ОЮ 0,005 О, ООО
0,044
0,019
0,011
10 20 Параметр останова
50
Рис. 19 - Функция Химмельблау: погрешность решения е , как функция
параметра останова т
Рисунки 18, 19 показывают, что для функции Химмельблау с увеличением параметра останова т погрешность найденного решения е уменьшается, одновременно возрастает число итераций te.
5.3. Функция Растригина (Rastrigm). Рассмотрена двумерная функция, обратная функции Растригина (рисунок 20)
F(.х, у) = -(20 + x2 + У2 - 10(cos2лx + cos2яy))
в области В = {(х,у)| х е[-2,2],у е[-2,2]}. Глобальный максимум этой функции достигается в точке (0; 0) и равен нулю. Если не указано иное, размер области локального поиска А полагается равным 0,08.
X
а) б)
Рис. 20 - Поверхность (а) и линии уровня (б) двумерной функция Растригина
Варьирование числа агентов-разведчиков. Исследование выполнено при числе агентов-разведчиков £, равном 1, 2, 4, 8, 16, 32. Результаты исследования иллюстрируют рисунки 21 - 23.
Рисунок 23 представляет оценку вероятности локализации глобального экстремума р1, которая определяется на основе 30 прогонов программы с разными начальными положениями агентов, как относительное число решений, находящихся в области «притяжения» глобального экстремума.
Рисунки 20 - 23 показывают, что для функции Растригина с ростом числа агентов-разведчиков £ число итераций te не имеет четкой тенденции, погрешность решения е достигает максимального значения при £ = 8,
вероятность локализации глобального экстремума р1, как и следовало ожидать, возрастает.
Рис. 21 - Функции Растригина: число итераций te в функции числа агентов-
разведчиков 5
Рис. 22 - Функции Растригина: погрешность решения е, как функция числа
агентов-разведчиков 5
Варьирование размера области локального поиска. Исследование выполнено при параметре Д , равном 0,04; 0,08; 0,2; 0,4. Результаты исследования представлены на рисунках 24 - 26.
Рис. 23 - Функции Растригина: оценка вероятности локализации глобального максимума р1 в функции числа агентов-разведчиков 5
0,04 0,08 0,2 0,4
Размер области локального поиска
Рис. 24 - Функции Растригина: число итерации te в функции размера
области локального поиска Д
Рис. 25 - Функции Растригина: погрешности решения е, как функция размера области локального поиска Д
100,0 96,7
0,04 0,08 0,2 0,4
Размер области локального поиска
Рис. 26 - Функции Растригина: оценка вероятности локализации глобального экстремума р1 в функции размера области локального поиска А
Рисунки 24, 25 показывают, что для функции Растригина с ростом размера области локального поиска А число итераций te медленно убывает, а погрешность решения е быстро возрастает. Из рисунка 26 следует, что вероятность локализации глобального экстремума р1 с ростом данного параметра достигает минимального значения при А = 0,08 и с дальнейшим ростом А существенно увеличивается.
Варьирование параметра останова т. Исследование выполнено при значениях указанного параметра, равных 5, 10, 20, 50. Результаты исследования иллюстрируют рисунки 27 - 29.
Рисунки 27 - 29 показывают, что для функции Растригина с ростом параметра останова т число итераций te быстро растет, погрешность полученного решения е убывает, а вероятность локализации глобального экстремума р1 медленно и почти линейно возрастает.
5.4. Сравнение результатов исследования. Приведенные выше результаты исследования объединены на рисунках 30 - 32, на которых ромбами отмечены результаты, относящиеся к функции Розенброка,
треугольниками Растригина.
- к функции Химмельблау, квадратами - к функции
Рис. 27 - Функции Растригина: число итераций te в функции параметра
останова т
Рис. 28 - Функции Растригина: погрешности решения s, как функция
параметра останова т
Рис. 29 - Функции Растригина: оценка вероятности локализации глобального экстремума р1 в функции параметра останова т
г 0.045 0.04 0.035
10
20
а)
30
40 Ш
10
20
б)
30
40 Ш
Рисунок 30 - Число итераций te (а) и погрешность решения е (б) в функции
числа агентов-разведчиков 5
Рисунок 31 - Число итераций te (а) и погрешность решения е (б) в функции
размера области локального поиска А
Рисунок 30 показывает, что с увеличением числа агентов-разведчиков 5 для одноэкстремальной функции Розенброка имеет место уменьшение числа итераций te, а для многоэкстремальных функций Химмельблау и Растригина - небольшое увеличение этого числа. Аналогично, погрешность решения е для функций Розенброка и Химмельблау уменьшается с ростом 5, а для функции Растригина - незначительно увеличивается. Асимптотически, для
всех трех указанных функций метод обеспечивает близкие результаты; начиная с числа агентов-разведчиков, равного £«30, наблюдается стагнация величин te, е .
Рисунок 32 - Число итераций te (а) и погрешность решения е (б) в функции
параметра останова т
Рисунок 31. позволяет сделать следующие выводы. Для функции Розеброка и Химмельблай с ростом размера области локального поимка А имеет место уменьшение числа итераций te. Погрешность решения е при
этом возрастает с существенно различной скоростью для разных функций.
Из рисунка 32 следует, что с ростом параметра останова т число итераций te растет для всех трех функций с примерно одинаковой скоростью. В то же время, погрешности е быстро убывают и при т ~ 40 становятся близкими.
Заключение
Выполненное исследование показало, что метод роя пчел обеспечивает 100 % вероятность локализации глобального экстремума функций Розенброка и Химмельблау; для функции Растригина, которая является многоэкстремальной и имеет сложную топологию, та же вероятность зависит от значений свободных параметров метода и изменяется от 20 % до 97 %.
В развитие работы планируется исследование эффективности метода в условиях высокой размерность вектора варьируемых параметров, а также более широкое исследование влияния свободных параметров метода на его эффективность. Планируется также разработка и исследование модификаций метода, ориентированных на параллельные вычислительные системы различной архитектуру (кластеры, системы с общей памяью, графические процессорные устройства). Наконец, имеется в виду сравнение эффективности и интеграция различных поведенческих методов.
Литература
1. Weise T. Global Optimization Algorithms - Theory and Application: Ph.D. thesis / University of Kassel.- 2008.
2. Beni G., Wang J. Swarm Intelligence // Annual Meeting of the Robotics Society: Proceedings of Seventh International Conference.- Tokyo: RSJ Press, 1989.- pp. 425-428.
3. Bonabeau E., Dorigo M., Theraulaz G. Swarm Intelligence: From Natural to Artificial Systems. - New York: Oxford University Press, 1999. - 320 p.
4. Karaboga D. An idea based on honey bee swarm for numerical optimization. Technical report - TR06, Erciyes University, Engineering Faculty, Computer Engineering Department 2005.
5. Pham D.T., Ghanbarzadeh A., Koc E, Otri S, Rahim S., Zaidi M. The Bees Algorithm. Technical Note, Manufacturing Engineering Centre, Cardiff University, UK, 2005
6. Pham D.T., Ghanbarzadeh A., Koc E., Otri S., Rahim S., Zaidi M. The Bees Algorithm - A Novel Tool for Complex Optimisation Problems Manufacturing Engineering Centre, Cardiff University, Cardiff CF24 3AA, UK.
7. Python Programming Language - Official Website [Электронный ресурс] // (http://pkolt.ru/pages/python/).
8. Tang K., Yao X., Suganthan P.N., MacNish C., Chen Y.P., Chen C.M., Yang Z. Benchmark Functions for the CEC'2008 Special Session and Competition on Large Scale Global Optimization. - Nature Inspired Computation and Applications Laboratory, USTC, China, 2007.
9. Chong S.C., Low M.Y.H. A Bee Colony Optimization Algorithm to Job Shop Scheduling // Winter Simulation Conference: Proceedings of the 38th Conference on Winter simulation. -Monterey: Monterey Press, 2006. - P. 1954-1961.
10. Lucic P., Teodorovic D. Bee System: Modeling Combinatorial Optimization Transportation Engineering Problems by Swarm Intelligence // Transportation Analysis: Proceedings of the Triennial Symposium TRISTAN IV. - Sao Miguel: Azores Press, 2001. - P. 441445.
11. Teodorovic D., Dell'Orco M. Bee Colony Optimization - a Cooperative Learning Approach to Complex Transportation Problems // Advanced OR and AI Methods in Transportation: Proceedings of 16th Mini-EURO Conference and 10th Meeting of EWGT (13-16 September 2005). -Poznan: Publishing House of the Polish Operational and System Research, 2005. - P. 51-60.
12. Nakrani S., Tovey C. On Honey Bees and Dynamic Allocation in an Internet Server Colony // Adaptive Behavior. - 2004. - №12. - P. 223240.
13. Quijano N., Passino K.M. Honey Bee Social Foraging Algorithms for Resource Allocation: Theory and Application. - Columbus: Publishing House of the Ohio State University, 2007. - 39 p.
14. Sumpter D.J.T., Broomhead D.S. Formalising the Link between Worker and Society in Honey Bee Colonies // Lecture Notes in Computer Science: Proceedings of the First International Workshop on Multi-Agent Systems and Agent-Based Simulation.- MABS '98 LNAI, 1998.-pp. 95-110.
15. Олейник Ал.А. , Субботин С.А. Интеллектуальный метод мультиагентного поиска поиска в многомерном пространстве с использованием прямой связи между агентами // Складш системи i процеси, 2008, №2.- с. 102 108.