НАУКА и ОБРАЗОВАНИЕ

УДК 519.6

Россия, МГТУ им. Н.Э. Баумана

apkarpenko@mail.ru

dmitry.shurov@mail.ru

Введение

Искусственные иммунные системы (ИИС) – это информационная технология, использующая понятия, аппарат и некоторые достижения теоретической иммунологии для решения прикладных задач. На основе ИИС разработано большое число различных методов для численного решения прикладных задач анализа данных, распознавания образов, компьютерной безопасности, диагностики неисправностей технических объектов, оптимизации и т.д.

Методы ИИС, ориентированные на решение задачи глобальной оптимизации, вдохновлены некоторыми аспектами поведения иммунной системы человека в процессе защиты ею организма от внешних факторов (патогенов и антигенов). Защитные клетки иммунной системы (антитела) при этом претерпевают множество изменений, целью которых является создание клеток, обеспечивающих наилучшую защиту от данного фактора. В результате создается большое число различных антител, и в борьбе с патогеном побеждают те из них, которые оказались наилучшим образом приспособлены для защиты. Иммунная система человека сохраняет в своей «памяти» победившие антитела, на основании чего именно такие антитела производятся при повторном проникновении в организм схожего патогена.

При разработке методов оптимизации на основе ИИС не ставится цель в точности воспроизвести детали функционирования живых иммунных систем. Авторы таких методов ограничиваются применением лишь некоторых принципов функционирования ИИС для конструирования механизмов оптимизации. Например, степень близости антител и антигенов формализуют в виде евклидовой нормы в пространстве варьируемых параметров, степень приспособленности антитела для борьбы с данным антигеном – в виде инвертированного значения целевой функции и т.д.

Основной целью публикации является представление оригинального гибридного метода глобальной оптимизации на основе ИИС. Метод получил наименование SIA (SubplexImmuneAlgorithm). В качестве базового метод использует известный метод HIA [1]. Модификация заключается, прежде всего, в отказе от процедуры мутации, как средства локального поиска. Вместо этого метод SIA использует для локального поиска алгоритм SUBPLEX  [2], благодаря чему гарантируется, что до завершения полной проверки из популяции не будут исключены изначально кажущиеся плохими решения. Научная новизна работы заключается в новизне предложенного метода глобальной оптимизации.

Работа организована следующим образом.

В первом разделе дана постановка задачи глобальной оптимизации и общая схема иммунных методов ее решения.

Второй раздел содержит обзор значительного числа методов глобальной оптимизации, построенных на основе ИИС, в частности, в этом разделе дано описание базового метода HIA.

Третий раздел представляет разработанный метод SIA. Здесь же на основе сравнительного исследования эффективности алгоритмов локальной оптимизации, реализованных в известной программной библиотеке NLopt [3], обосновано использование алгоритма SUBPLEX.

Четвертый раздел посвящен разработке программного обеспечения, реализующего метод SIA. Реализация выполнена на языке программирования общего назначения C++. Разработанное программное обеспечение SIASearchпредставляет собой многофункциональную расширяемую систему для тестирования и исследования эффективности методов оптимизации. В этом же разделе с использованием известных тестовых функций Розенброка, Растригина, и Химмельблау [4] выполнено тестирование SIASearch, показавшее корректность его функционирования.

В пятом разделе представлены результаты исследования эффективности алгоритма SUBPLEX, метода SIA, а также результаты сравнительного исследования эффективности методов HIA, SIA и методов глобальной оптимизации CRS, MLSL, ISRES из библиотеки NLopt.

В заключении сформулированы основные результаты работы и очерчены перспективы ее развития.

1. Постановка задачи и общая схема методов иммунной оптимизации

Рассматриваем задачу глобальной условной минимизации скалярной целевой функции  вещественных переменных  в области допустимых значений :

.                                         (1)

Здесь  ‑ размерность вектора . Общее число ограничений, определяющих область , обозначаем . Если далее не оговорено противное, то имеется в виду частный случай задачи (1) ‑ задача безусловной оптимизации, когда . Полагаем, что для этой задачи задана область поиска в виде параллелепипеда

,

где ,  ‑  заданные константы.

В работе используем следующие основные понятия:

·  клетка, антитело ― решение ;

·  аффинность ― величина, обозначающая полезность клетки и равная соответствующему значению целевой функции, взятому с противоположным знаком;

·  популяция― множество клеток;

·  клон ― клетка, образованная из другой клетки путем ее полного копирования;

·  мутация ― случайное изменение компонентов вектора .

Используем также следующие обозначения указанных и других сущностей:

·  ― клетка;

·   ― значение целевой функции, соответствующее данной клетке;

·  ― популяция клеток;

·  ― компоненты клетки ;

·  ― аффинность клетки;

·  ― популяция клеток памяти (не являющихся частью популяции );

·   ― множество клонов клетки , созданных на текущей итерации метода;

·   ― множество клонов всех клеток текущей популяции;

·   ― ‑й клон клетки ;

·   ― лучший клон клетки , то есть клон, которому соответствует минимальное значение целевой функции;

·   ― число клонов, порождаемых каждой из клеток популяции ;

·   ― минимальное значение целевой функции, найденное методом численной оптимизации;

·   ― среднее по мультистартам значение целевой функции, найденное этим методом;

·   ― среднее по мультистартам число вычислений целевой функции (число испытаний).

Во введенных обозначениях основные особенности методов иммунной оптимизации можно сформулировать следующим образом.

1) Все методы иммунной оптимизации работают с популяцией клеток

,

в которой клетка  представляет собой вектор в пространстве варьируемых переменных.

2) Над клетками популяции осуществляют процедуры клонирования и мутации. При клонировании клетки  создают множество клеток , идентичных исходной:

.                                  (2)

При мутации компоненты вектора  изменяют случайным образом по схеме

,                                      (3)

где  ― шаг мутации, способ вычисления которого различен для различных методов иммунной оптимизации.

3) В результате процедуры клонального отбора лучшие из потомков замещают родителей по формуле

                                     (4)

4) Процедура сжатия популяции устраняет избыточность, исключая худшее из решений , удовлетворяющих условию

,                                             (5)

где  ‑ символ некоторой меры близости клеток, например, евклидова норма;  ‑ порог сжатия.

5) Множество клеток памяти  хранит некоторое число наилучших решений, найденных за все время работы метода.

Процедуры клонирования и мутации обеспечивают локальную оптимизацию в окрестности каждой из клеток. Процедура сжатия популяции преследует цель ускорения сходимости метода и сокращения вычислительных затрат.

2. Обзор методов иммунной оптимизации

          Последовательно рассматриваем методы CLONALG, opt-AiNet, BCA, HIA, I-opt-AiNet и, наконец, метод T-Cell Model.

2.1. Метод CLONALG

Метод CLONALG [5] был разработан для решения задач машинного обучения и распознавания образов. Позднее метод адаптировали для задач оптимизации. Основные принципы метода: наличие клеток памяти, отбор и клонирование наиболее полезных антител, исключение наименее полезных антител, выбор клонов пропорционально их полезности, обеспечение и сохранение разнообразия популяции антител.

Общая схема метода CLONALG имеет следующий вид [5].

1) Определяем полезность  каждого из антител  множества ; .

2) Во множестве  оставляем некоторое число антител с наибольшей полезностью, а остальные антитела исключаем.

3) Клонируем антитела из модифицированного множества  пропорционально их полезности и помещаем полученные клоны во множество .

4) Множество клонов  подвергаем мутации, интенсивность  которой уменьшается с ростом полезности по формуле

,

где  ― положительный свободный параметр.

5) Определяем полезность всех мутировавших клонов множества , и находим на этой основе клетку с наивысшей аффинностью. Если полезность данной клетки выше, чем у соответствующей ему клетки памяти, то заменяем последнюю указанным мутировавшим клоном.

6) Осуществляем замену -ой части антител из множества  новыми клетками. Здесь  ‑ свободный параметр метода.

2.2. Метод opt-AiNet

Главной особенностью метода opt-AiNet является сохранение всех найденных локальных и глобальных оптимумов целевой функции. Метод использует динамический размер популяции, управляемый механизмом сжатия [6].

Схему метода opt-AiNet передает следующая последовательность его основных шагов.

1)  Инициализируем популяцию клеток.

2)  Определяем полезность  всех клеток популяции и нормализуем значениях этих величин так, чтобы нормализованная полезность  принадлежала интервалу ; .

3)  Порождаем  клонов для каждой из клеток.

4)  Осуществляем мутацию каждого клона пропорционально полезности его родительской клетки по формуле

,     ,

где  ― решение после мутации клетки ;  ― свободный положительный параметр метода.

5)  Определяем полезность мутировавших клонов.

6)  Для каждой родительской клетки и ее клонов выбираем наилучшее решение и исключаем из популяции остальные решения.

7) Вычисляем среднюю полезность популяции, и если она уменьшилась по сравнению с предыдущим шагом, то возвращаемся к шагу 2.

8) Определяем расстояние между всеми клетками популяции и осуществляем сжатие тех клеток, расстояние между которыми ниже заданного порога  (см. формулу (5)). Оставшиеся клетки помещаем во множество клеток памяти.

9) Вводим в популяцию -ю часть случайно сгенерированных клеток.

Здесь и далее  ‑ нормальная случайная величина cматематическим ожиданием a и среднеквадратическим отклонением b.

В работе [13] рекомендованы следующие значения свободных параметров метода: ; ;; . Начальный размер популяции  рекомендован равным 20. Гибридизация метода opt-AiNet с методом локального поиска рассмотрена в работе [7]. В этом случае к рассмотренной схеме метода добавляется этап, предполагающий локальную оптимизацию наилучшего решения.

2.3. МетодBCA (B-cell algorithm)

Метод BCA, как и метод opt-AiNet, основан на принципе клонального отбора [8]. Метод ориентирован на поиск глобального оптимума целевой функции, имеющий сложный ландшафт, и может обеспечить высокую эффективность поиска при небольших размерах начальной популяции . Особенностью метода является использование оператора так называемой смежной гипермутации (contiguoussomatichypermutation), который предполагает мутацию не отдельных случайных компонентов клетки, а целой области из нескольких смежных (соседних по номерам) компонентов.

Схема метода BCAпредставлена ниже [8].

1) Инициализируем популяцию клеток.

2) Вычисляем полезность всех клеток популяции.

3) Для каждой клетки  порождаем набор клонов , число клеток в котором равно начальному размеру популяции, то есть ; .

4) Из каждого набора клонов  случайно выбираем одну клетку и осуществляем случайное изменение ее компонентов аналогично методу CLONALG.

5) Ко всем клонам применяем оператор гипермутации.

5) Вычисляем полезность всех мутировавших клонов. Для каждого набора клонов определяем клетку с наибольшей полезностью и, если ее полезность превышает полезность родительской клетки, заменяем последнюю первой.

6) Если условие окончания итераций не выполнено, то переходим к шагу 2.

2.4. Метод HIA

Метод HIA (Hybrid Immune Algorithm) [1] можно рассматривать как модификацию метода opt-AiNet. В отличие от последнего метода, размер популяции метода HIA остается постоянным в процессе всего итерационного процесса; каждая из клеток имеет дополнительный атрибут, называемый возрастом клетки; мутацию клонов осуществляют двумя способами ‑ независимо для каждой компоненты клетки, либо только для одной ее компоненты. Метод HIA применяет стратегию, в соответствие с которой новые клетки создаются не на основе равномерного распределения по всей области поиска, как в методе opt-AiNet, а с использованием распределения Гаусса с динамически изменяющимися параметрами. Это приводит к постепенному сужению области поиска, но, благодаря второму типу мутации, не препятствует широкому исследованию пространства поиска.

В соответствии с работой [1] схему алгоритма HIA можно представить в следующем виде.

1) Инициализируем популяцию клеток.

2) Осуществляем клонирование, при котором каждая клетка популяции порождает одинаковое число  клонов.

3) Подвергаем клоны мутации.

4) Для каждой родительской клетки определяем лучший клон, и если полезность этого клона превосходит полезность родителя, заменяем родительскую клетку этим клоном. В противном случае возраст родителя увеличиваем на единицу.

5) Клетки, возраст которых достиг заданного значения перемещаем из основной популяции во множество клеток памяти (п.3).

6) Осуществляем сжатие клеток памяти и, если число этих клеток после сжатия превышает заданную величину, исключаем худшие из них.

7) Если не выполнен критерий окончания итераций, то переходим к шагу 2.

Первый способ мутации определяет формула

,

где

Второй способ мутации осуществляется в 20 % случаев и заключается в замене одной случайно выбранной компоненты  родительской клетки  на случайное число, равномерно распределенное в заданном интервале ; , .

2.5. Метод I-opt-AiNet

Метод I-opt-AiNet представляет собой модификацию метода opt-AiNet, которая преследует цель уменьшить число испытаний в процессе поиска решения [9]. Изменения касаются определения числа выбираемых для клонирования решений и способа их выбора, числа удаляемых из популяции и добавляемых в нее решений, способа выбора решений для клонирования и числа порождаемых ими клонов, а также критерия остановки. Метод использует понятие лучшей средней полезности по всем предыдущим итерациям (Best Fitness Average, BFA), обозначаемой . Для того, чтобы предотвратить повторные ресурсоемкие вычисления полезности клеток, в процессе работы метода те клетки, значение полезности которых однажды уже было вычислено, помечают как «вычислены».

Схема метода I-opt-AiNet имеет следующий вид [9].

1) Инициализируем популяцию клеток.

2) Вычисляем полезность всех клеток популяции и эти клетки помечаем как «вычислены».

3) Из популяции исключаем клетки, полезность которых ниже величины .

4) Осуществляем клонирование клеток пропорционально их полезности, так что

,    ;    .

Здесь  ― максимальное число клонов, которое может быть порождено клеткой; , ― максимальное и минимальное текущие значения полезности среди всех клеток популяции соответственно;  ― нормированная полезность клонируемой клетки.

5) Осуществляем мутацию клонов по схеме метода opt‑AiNet. Мутацию принимаем лишь в том случае, если результирующая клетка расположена в рассматриваемой части области поиска.

6) Вычисляем полезность всех клонированных клеток и помечаем их как «вычислены».

7) Из полученной популяции исключаем те клетки, полезность которых ниже .

8) В каждой паре клеток, расстояние между которыми меньше заданной величины , исключаем клетку с меньшей полезностью.

9) В популяцию добавляем случайно сгенерированные клетки, число которых равно числу клеток, исключенных на шаге 3.

10) Если не выполнен критерий окончания итераций, то переходим к шагу 2.

В работе [9] рекомендованы следующие значения свободных параметров метода: максимальное число клонов каждой клетки ; порог сжатия ; начальное число клеток .

2.6. Метод T-Cell Model

В отличие от рассмотренных выше методов, метод T-Cell Model предназначен для решения задач глобальной условной оптимизации [10]. Метод использует три субпопуляции клеток ‑ свободные клетки , исполнительные клетки  и клетки памяти M. Свободные клетки обеспечивают разнообразие поиска и представляются двоичными строками в коде Грея. С помощью исполнительных клеток, также представляемых в коде Грея, исследуют части пространства поиска, лежащие вне области поиска, но недалеко от ее границы. Клетки памяти кодируют вещественными числами и используют для организации локального поиска в окрестностях лучших найденных решений.

Для того чтобы повысить вероятность локализации глобального экстремума целевой функции, метод T-Cell Model  использует так называемый динамический допуск  (Dynamic Tolerance Factor, DTF). Величину  назначают для каждой из трех указанных выше субпопуляций по следующему правилу:

‑ вычисляем для каждой из субпопуляций величину , равную сумме значений ограничивающих функций, для которых клетками данной субпопуляции нарушено хотя бы одно из ограничений;

‑ полагаем ,  для субпопуляции ,  соответственно, где  ‑ числа клеток в этих популяциях;

‑  если для субпопуляции  имеет место неравенство , то принимаем для этой субпопуляции ;

‑ в той же ситуации для субпопуляции  принимаем ;

‑ для субпопуляции Mво всех случаях принимаем .

Для субпопуляции  мутацию клеток не используют. Для клеток субпопуляций , M, используют следующие законы мутации.

1.1) Вариант , . Определяем величину наиболее нарушенного ограничения, и если она превышает , то каждый бит клетки pинвертируем с вероятностью . В противном случае по тому же правилу инвертируем биты только одной случайно выбранной компоненты клетки.

1.2) Вариант , . Каждый бит клетки pинвертируем с вероятностью .

2) Вариант . Мутацию клетки pпроизводим по формуле

,

где  ― случайное число, равномерно распределенное в интервале ; ― целое число (свободный параметр метода).

Важной в методе T-Cell Model  является процедура отсеивания клеток, обеспечивающая сохранение лучших найденных решений. Для клеток, находящихся в области допустимых значений D, при отсеивании учитывают значения целевой функции, а для клеток вне этой области ― величину . Любая из клеток, принадлежащих области D, считается лучше любой из клеток вне этой области.

Схема метода T-Cell Model  состоит из следующих основных шагов.

1) Случайно генерируем заданное число свободных клеток.

2) Вычисляем значение динамического допуска  для популяции свободных клеток.

3) С учетом динамического допуска определяем, какие из свободных клеток принадлежат области D, а какие на принадлежат.

4) Заданный процент свободных клеток включаем в популяцию исполнительных клеток.

5) Вычисляем значение динамического допуска для популяции исполнительных клеток.

6) Повторяем 50 раз процедуру мутации исполнительных клеток и определяем, какие из исполнительных клеток лежат в области Dс учетом их динамического допуска, а какие не лежат.

7) Заданный процент исполнительных клеток включаем в популяцию клеток памяти.

8) Повторяем 100 раз для клеток памяти действия, указанные в 6.

9) Описанные выше шаги выполняем заданное число раз.

2.7. Резюме

Во всех рассмотренных методах иммунной оптимизации задача интенсификации поиска (задача поиска локальных минимумов) решается с помощью операторов клонирования и мутации, которые на основе уже существующих в иммунной системе клеток продуцируют множество дочерних клеток. Задача диверсификации поиска (исследования пространства поиска) решается, преимущественно, путем введения в область поиска новых случайно сгенерированных решений. С целью сокращения вычислительных затрат, некоторые из рассмотренных методов имеют средства устранения близких решений. Наконец, все рассмотренные методы используют множество «клеток памяти», содержащее лучшие найденные решения.

Сравнительная характеристика рассмотренных методов представлена в таблице 1.

Вероятно, самой сильной стороной методов оптимизации на основе ИИС является наличие механизма памяти. Клетки памяти хранят информацию о найденных ранее наилучших решениях и позволяют использовать ее для повышения эффективности глобального поиска. Однако, на наш взгляд, рассмотренные методы не в полной мере используют потенциальные преимущества этого механизма, поскольку новые клетки в этих методах включаются в популяцию случайно по всей области поиска, независимо от состояния памяти. Единственным исключением в этом смысле является метод HIA, в котором с увеличением номера итерации информация о ландшафте целевой функции, хранящаяся в клетках памяти, все в большей степени влияет на генерацию новых решений. Благодаря этому неперспективные области пространства параметров исключаются из рассмотрения, и поиск постепенно концентрируется в окрестности наилучших из найденных решений, что позволяет значительно повысить эффективность глобального поиска.

Таблица 1 ― Сравнение методов иммунной оптимизации

Другой важной отличительной особенностью рассмотренных методов является использование ими процедуры сжатия множества решений, которая устраняет клетки, расположенные в пространстве поиска слишком близко к другим клеткам. Такая процедура позволяет уменьшить избыточность популяции и, в результате, повысить скорость поиска. Поскольку размер популяции обычно невелик, вычислительные затраты на сжатие оказываются небольшими и оправдываются за счет уменьшения числа испытаний.

Главная слабость рассмотренных методов оптимизации состоит в низкой эффективности используемой ими процедуры локального поиска на основе мутации клеток. Эта процедура, по сути, представляет собой один из вариантов случайного поиска, эффективность которого, как хорошо известно, резко снижается с увеличением размерности задачи. Если даже начальное приближение сгенерировано в области притяжения глобального минимума, данная процедура может привести к потере этого решения, поскольку другие решения, находящиеся в областях притяжения локальных минимумов, могут иметь лучшие значения целевой функции.

3. Метод SIA

Разработка метода SIA имела целью преодолеть главный недостаток рассмотренных в предыдущем разделе методов иммунной оптимизации (низкую эффективность процедуры локального поиска) путем использования современного высокоэффективного алгоритма локальной оптимизации. Кроме того, метод SIA использует оригинальную адаптивную процедуру генерации начальных приближений, которая учитывает информацию о ранее найденных локальных решениях, хранящихся в клетках памяти.

Для выбора алгоритма локального поиска было проведено сравнительное тестирование алгоритмов, представленных в программной библиотеке NLopt [3]. В качестве тестовой функции использована овражная функция Розенброка (п.4), имеющая в точке минимума значение . При тестировании был применен один из стандартных для библиотеки NLopt критерий окончания итераций:  для всех . Здесь  ‑ текущий шаг по i-ой компоненте вектора варьируемых параметров;  ‑ текущее значение этой компоненты. Использовано число мультистартов, равное 100. Рассмотрены три критерия эффективности алгоритма:

‑ минимальное достигнутое значение целевой функции ;

‑ среднее по мультистартам значение указанного критерия ;

‑ среднее по мультистартам число испытаний .

Из результатов тестирования, представленых в таблицах 2 – 4, следует, что алгоритмы NEWUOA, SUBPLEX показывают примерно одинаково хорошие результаты. Однако поскольку алгоритм NEWUOAможет иметь низкую эффективность в случае, когда вторая производная целевой функции не существует [3], для гибридизации нами был выбран алгоритм SUBPLEX.

Таблица 2 ― Сравнение методов локального поиска по критерию

Таблица 3 ― Сравнение методов локального поиска по критерию

Алгоритм SUBPLEX (subspace‑searching simplex) предназначен для оптимизации сложных зашумленных недифференцируемых функций с большим числом варьируемых параметров. Метод разработан на основе симплекс-метода Нелдера-Мида [11] и использует его сильные стороны, устраняя недостатки. SUBPLEX осуществляет разбиение вектора варьируемых параметров xна подвекторы и применяет алгоритм Нелдера-Мида независимо к каждому из этих подвекторов. Таким образом, метод SUBPLEX осуществляет поиск в совокупности подпространств пространства варьируемых параметров  (см. ниже).

Таблица 4 ― Сравнение методов локального поиска по критерию

Схема метода SIA имеет следующий вид.

1) Инициализируем популяцию с заданным числом клеток.

2) Методом SUBPLEX для каждой из клеток популяции выполняем локальный поиск и заменяем координаты этой клетки на координаты, полученные в результате локального поиска.

3) Перемещаем все клетки популяции во множество клеток памяти, исключаем их из популяции, и добавляем в популяцию новые клетки с координатами

,

где

;

;

.

Здесь  ‑ значение коэффициента обучения на текущей итерации t;  ‑ клетки текущего множества памяти;  и  ― коэффициенты, задающие скорость обучения, рекомендованные значения которых равны ;  [1].

4) Осуществляем сжатие множества клеток памяти: для каждой пары этих клеток вычисляем расстояние между ними; если это расстояние меньше заданного порога , то исключаем клетку с худшим значением целевой функции; если число клеток памяти после сжатия превышает величину , то оставляем  лучших из них.

5) Если условия окончания итераций выполнены, то завершаем вычисления. В противном случае переходим к шагу 2.

В завершении раздела рассмотрим метод SUBPLEX. Метод использует следующие свободные параметры:

 ‑ коэффициент отражения ();

 ‑ коэффициент расширения ();

 ‑ коэффициент сужения ();

 ‑ коэффициент сжатия ();

 ‑ коэффициент уменьшения симплекса ();

 ‑ коэффициент уменьшения шага ();

 ― начальное значение компонентов вектора шагов;

,  ‑ минимальная и максимальная размерности подпространств, подчиненные ограничениям

,    ,

где  ‑ ближайшее целое большее .

В работе [2] рекомендованы следующие значения указанных параметров: ; ; ; ; ; ; ; .

Последовательность основных шагов метода SUBPLEX имеет следующий вид.

1) Определяем величины компонентов вектора шагов.

2) Выделяем подпространства в пространстве поиска.

3) Для каждого из подпространств применяем классический метод поиска Нелдера-Мида с так называемым внутренним критерием окончания итераций.

4) Проверяем выполнение внешнего условия окончания итераций. Если это условие не выполнено, то переходим к шагу 1. В противном случае завершаем вычисления.

Определение величины компонент вектора шагов. Компоненты этого вектора определяют величину и ориентацию начальных симплексов. На первой итерации все эти компоненты инициализируем значением . На последующих итерациях сначала масштабируем вектор шагов  по формуле

где  ― вектор приращения решения по сравнению с предыдущей итерацией;  ― число подпространств, выделенных на предыдущей итерации; . Затем корректируем компоненты вектора шагов по формуле

Выделение подпространств в пространстве поиска. В результате выполнения этого этапа определяем число подпространств, размерность каждого из подпространств, а также координатные оси, задающие эти подпространства. Основная идея формирования подпространств состоит в том, чтобы направление вектора шагов находилось преимущественно в первом подпространстве. Число и размерности подпространств должны удовлетворять условиям

,    ,

где  ― размерность -го подпространства.

Первым этапом определения подпространств является сортировка компонентов вектора шагов по убыванию их величин. Обозначаем исходный и отсортированный векторы шагов ,  соответственно, где .

На втором этапе определяем размерность  первого подпространства как целое число , которое обеспечивает максимальное значение функции

при условии

.

Выполнение этого условия гарантирует, что на основе оставшейся части вектора xтакже могут быть сформированы следующие подпространства.

По аналогичной схеме определяем ,  и т.д.

Поиск в подпространстве выполняем с помощью классического метода Нелдера-Мида до тех пор, пока размер симплекса не уменьшится, по меньшей мере, в  раз по сравнению с его начальным размером (внутренний критерий останова). В качестве размера симплекса используем евклидово расстояние между наилучшей и наихудшей его вершинами.

Проверка внешнего условия окончания итераций. Локальный поиск по методу SUBPLEXостанавливаем, когда выполнен внешний критерий остановки

,

где  ― заданная погрешность; .

4. Программная реализация

При выборе программных средств для реализации метода SIA были рассмотрены следующие варианты:

‑ реализация на языке программирования общего назначения;

‑ реализация на специализированном языке одного из математических пакетов.

Последний вариант упрощает и ускоряет разработку вычислительных программ, но имеет несколько недостатков, главные из которых – это невозможность осуществления низкоуровневой оптимизации программы, зависимость от среды выполнения, возможное снижение скорости вычислений при интерпретации программы средой выполнения. Реализация на стандартном языке программирования общего назначения позволяет обеспечить максимальную производительность вычислений за счет оптимизации программы, не требует установки дорогостоящих программных пакетов и, что очень важно, обеспечивает возможность включения разработанной программы в состав другого программного обеспечения в виде библиотечной процедуры. На основе этих соображений было принято решение о реализации метода SIA на языке программирования C++, как одном из наиболее распространенных языков для решения подобных задач.

Разработанная программа, названная нами SIA Search, представляет собой расширяемую систему для тестирования и исследования эффективности методов оптимизации. Основными функциями SIA Search являются:

‑ задание множества значений свободных параметров метода, интересующих исследователя;

‑ автоматическое исследование эффективности метода в режиме мультистарта при всех заданных значениях свободных параметров;

‑ возможность введения в программу новой целевой функции без перекомпиляции программы;

‑ возможность расширения программы новыми методами оптимизации.

Поясним используемые в программе SIA Search правила, на основании которых задается указанное выше множество значений свободных параметров метода.

Пусть  ‑ вектор интересующих исследователя свободных параметров рассматриваемого метода оптимизации, а

  ‑

‑ параллелепипед допустимых значений компонентов этого вектора.

Множество подлежащих исследованию значений свободных параметров должно быть задано в виде

,

где  ‑ набор допустимых значений свободного параметра ;  ‑ число значений этого параметра; .

Программа SIA Search для каждой из возможных комбинаций  выполняет мультистарт с заданным числом запусков, для каждого запуска определяет значения указанных в п.3 критериев оптимальности метода, по завершении мультистарта находит средние значения этих критериев и представляет результаты исследования в виде таблицы.

Введение в программу SIA Search новой целевой функции осуществляется с помощью динамической библиотеки. Единственное условие для создания библиотеки ‑ это наличие в операционной системе компилятора GCC.

4.1. Структура классов программы

Программа состоит из нескольких классов, структуру которых иллюстрирует рисунок 1. Сплошной линией на рисунке обозначено отношение наследования, а пунктиром ‑ отношение «использует».

1) Класс OptFunction содержит информацию о целевой функции. Его конструктор принимает строковое имя целевой функции в виде экземпляра вспомогательного класса Name, размерность целевой функции и указатель на функцию языка программирования Cи, тип которой определен с помощью оператора typedefв виде

typedefdouble (*f)(unsignedintn, double* x).

Рисунок 1 ‑ Структура основных классов программы SIASearch

Класс OptFunction включает в себя следующие общедоступные функции:

·  doublevalue(vector<double> x) ‑ вычисляет значение целевой функции в точке, заданной вектором x;

·  unsignedintdimension() ‑ возвращает размерность целевой функции;

·  Namename() ‑ возвращает имя целевой функции.

2) Класс OptParams осуществляет смену значений параметров, для которых осуществляется исследование эффективности метода оптимизации. Каждый параметр исследуемого метода представлен набором значений, заданных классомOptParamVector. Интерфейс класса OptParamVector составляют конструктор, осуществляющий функции добавления параметров, и функция step (), реализующая последовательный перебор заданных значений параметров. Задание множества параметров осуществляется с помощью функции

doubleget(Namename, vector<double> values).

Получение текущего значения параметра осуществляется обращением по его имени с помощью функции

double get(const char* name).

3) Класс OptResult предназначен для формирования результатов работы метода оптимизации и включает в себя следующие функции:

· vector<double> coords() ‑ возвращает координаты наилучшего найденного решения;

· doublevalue() ‑ возвращает значение целевой функции, соответствующее наилучшему найденному решению;

· unsignedintiterations() ‑ возвращает число итераций метода, затраченных на поиск полученного решения;

· unsignedintevaluations() ‑ возвращает общее число вычислений целевой функции.

4) Виртуальный класс OptMethod является базовым классом для всех методов оптимизации. Его конструктор принимает класс целевой функции OptFunction и экземпляр класса OptParams, содержащий параметры используемого метода. Класс OptMethod имеет единственную общедоступную функцию start(), которая запускает метод оптимизации и по его завершении возвращает экземпляр класса OptResult с результатами оптимизации.

Класс OptMethod имеет защищенные функции, обязательные для определения в производных классах:

·  void setParams(OptParams p);

·  void init();

·  void iterate();

·  bool stop();

·  OptResult bestResult().

Функция setParams() предназначена для получения параметров метода по их строковым идентификаторам и инициализации внутренних переменных класса, соответствующих этим параметрам.

Функция init() предназначена для задания начальных условий, требуемых для работы метода. Данная функция вызывается один раз для данного набора параметров.

Функция iterate() осуществляет выполнение основной итерации рассматриваемого метода оптимизации. Реализация данной функции полностью зависит от этого метода.

Функция stop() проверяет критерий окончания итераций (стагнация в течение заданного числа итераций либо достижение заданного числа вычислений целевой функции).

Функция bestResult() возвращает лучший найденный на текущей итерации результат.

5) Класс Multistart предназначен для сбора, обобщения и вывода результатов исследования. В его задачи входит подсчет числа проведенных запусков исследуемого метода оптимизации, хранение результатов оптимизации по каждому запуску, осреднение этих результатов по всем проведенным стартам, вывод в файл в виде таблицы осредненных результатов и лучшего результата по всему мультистарту.

6) Каждая клетка искусственной иммунной системы является экземпляром класса Cell, который имеет функции, осуществляющие клонирование и мутацию клеток, вычисление расстояния между клетками, а также получение вектора координат клетки.

7) Класс AIS включает в себя функции, реализующие процедуры, специфические для методов иммунной оптимизации:

·  Cell randomCell() ‑ создает клетку, координаты которой генерируются случайно в пределах области поиска;

·  vector<Cell> clone(Cell& A, unsignedintNc) – порождает множество клонов заданной клетки и возвращает массив, содержащий эти клоны;

·  boolstagnation() ‑ определяет наличие ситуации стагнации, когда отсутствует улучшение значения целевой функции для лучшей из клеток памяти в течение заданного числа итераций;

·  OptResultbestResult() ― возвращает экземпляр класса OptResult, содержащий информацию о наилучшей из клеток памяти;

·  voidreplace(unsignedinti, constCell& cell) ― удаляет клетку  популяции и помещает на ее место клетку, заданную параметром cell;

·  voidselectBest(vector<Cell>& P, unsignedintN) – осуществляет сортировку популяции клеток по возрастанию значения целевой функции и сокращает популяцию до указанного размера;

·  voidsuppress(vector<Cell>& P, doubleD) – осуществляет процедуру сжатия популяции с порогом сжатия .

8) Класс SIA является основным классом, реализующим метод SIA. В данном классе определена функция iterate(), выполняющая все шаги метода; функция generateNew(), создающая новые клетки популяции на основе состояния клеток памяти; вспомогательная функция maxDeltaCoords(unsignedintcoord), определяющая максимальную разность координат клеток памяти.

9) Класс HIA был разработан для реализации сравнительного исследования эффективности метода SIA. Функции класса HIA близки к рассмотренным функциям класса SIA. Наиболее существенное отличие заключается в наличии функции mutate(vector<Cell>& P, Cell& p), которая осуществляет локальный поиск с помощью указанных в п.2.4 видов мутации.

4.2. Тестирование программы

Тестирование метода и программы SIA Search выполнено на функциях Розенброка, Химмельблау и Растригина [4].

Функция Розенброка является унимодальной, но овражной, и имеет вид

.

Область поиска ограничиваем гиперкубом , размеры которого равны . В этой области минимум  функция достигает в точке .

Функция Химмельблау

двумерна и имеет четыре локальных минимума, являющихся одновременно глобальными минимумами:

;

;

;

.

Важно, что все точки минимума функции, кроме первой, имеют иррациональные значения координат.

Функция Растригина по каждому из измерений представляет собой сумму квадратичной и косинусоидальной составляющих, благодаря чему функция имеет большое число регулярно расположенных локальных минимумов. Функцию определяет формула

.

Область поиска ограничивает гиперкуб, границы которого задают константы . Глобальный минимум  функция достигает в точке .

          Далее нам понадобятся также тестовые функции Швефеля и Шекеля.

          Функция Швефеля [4] является многоэкстремальной и обладает так называемым ложным глобальным минимумом. Многие методы поиска в этой ситуации склонны сходиться к последнему. Функцию Швефеля определяет формула

.

Поиск осуществляем в гиперкубе , в котором . Глобальный минимум функции в этом случае равен  и достигается в точке .

          Функция Шекеля [4] замечательна тем, что положение и глубину всех минимумов этой функции задает сам пользователь. Функция имеет вид

,

где  ‑ число минимумов;  ‑ вектор, определяющий координаты i–го минимума;  ‑ константа, задающая значение этого минимума.

Результаты тестирования программы SIA Search для двумерных вариантов функций Розенброка, Химмельблау и Растригина представлены в таблице 5. В верхних подстроках таблицы приведены значения, полученные с помощью метода SIA, а в нижних подстроках ‑ аналитические значения. Таблица показывает корректность метода и его программной реализации.

Таблица 5 — Результаты тестирования метода SIA:

5. Исследование эффективности

В данном разделе представлены некоторые результаты исследования эффективности алгоритма SUBPLEX, метода SIA, а также результаты сравнительного исследования эффективности методов HIA, SIA и алгоритмов глобальной оптимизации CRS, MLSL, ISRES, представленных в программной библиотеке NLopt.

5.1. Алгоритм SUBPLEX.

Результаты исследования эффективности алгоритма SUBPLEX представлены в таблице 6. Исследование выполнено для функции Розенброка, имеющей размерность от 50 до 800. Во всех случаях использован мультистарт из 30 запусков. Таблица показывает, что для всех размерностей каждый из запусков обеспечил высокую точность локализации минимального значения функции. Очень важно, что при этом наблюдается почти линейная зависимость числа испытаний  от числа варьируемых параметров .

Таблица 6 — Исследование эффективности алгоритма SUBPLEX:

функция Розенброка

5.2. Метод SIA.

Прежде было выполнено исследование эффективности метода при поиске глобального минимума функции Растригина размерностью 12, 16, 20, 24, 28, а также смещенной функции Растригина тех же размерностей. Область поиска, как и ранее, была ограничена гиперкубом , где . Глобальный минимум последней функции в этой области расположен в точке

.

Во всех случаях использован мультистарт из 100 запусков.

Результаты исследования представлены в таблице 7, которая показывает, что метод успешно находит точный глобальный минимум функции Растригина, имеющей размерность, не превышающую 20. При более высоких размерностях, когда число локальных минимумов становится чрезвычайно большим, метод находит локальные минимумы, близкие к точке глобального минимума.

Таблица 7 — Результаты тестирования метода SIA:

функция Растригина;

Параметрическое исследование метода SIA выполнено для функций Растригина и Шекеля, имеющих размерности 16, 24, 30, 36, 42. В качестве варьируемых параметров рассмотрены порог сжатия , максимальное число клеток памяти , коэффициент обучения , коэффициент . Использован фиксированный размер популяции . В каждом исследовании варьировался только один из рассматриваемых свободных параметров, а остальные параметры принимали свои «базовые» значения ; ; ;  для функции Шекеля и  для функции Растригина. Число мультистартов принято равным 100. В качестве критерия окончания итераций во всех случаях использована стагнация в течение 20 итераций.

Результаты исследования показали, что только параметр  оказывает существенное влияние на эффективность метода. На этом основании приводим зависимости рассматриваемых критериев оптимальности только от этого параметра (рисунок 2).

Результаты исследования показывают, что среднее по мультистарту значение целевой функции  уменьшается (улучшается) с ростом числа клеток памяти  до величины, примерно равной 90, которая слабо зависит от размерности этой функции (рисунок 2б). С другой стороны, с увеличением  увеличивается среднее число испытаний . Таким образом, как и следовало ожидать, увеличение значений параметра  повышает качество поиска ценой увеличения числа испытаний.

а)                                                      б)

в)

► ‑ , ▲ ‑ , ▼ ‑ , ♦ ‑ , ■ ‑

Рисунок 2 ‑ Эффективность метода SIA:

функция Растригина; варьируемый параметр

Для функции Шекеля зависимость эффективности метода SIA от числа клеток памяти имела иной характер. Наилучшие результаты для всех указанных выше размерностей целевой функции были достигнуты в этом случае при , а для всех остальных  были получены примерно одинаковые результаты.

Разный характер зависимости критериев эффективности от параметра  для функции Растригина и Шекеля объясняется существенным различием ландшафтов этих функций.

5.4. Методы SIA, HIA, CRS, MLSL, ISRES

Сравнительное исследование эффективности указанных методов выполнено для функций Растригина и Шекеля, имеющих размерность 16, 24, 32, 40, 50. Для функции Растригина использована область поиска, указанная в п. 5.3. В случае функции Шекеля число минимумов  принято равным размерности целевой функции (); положения минимумов определяют константы

, ;

глубину минимумов задают величины , . Область поиска для функции Шекеля ограничена гиперкубом со стороной, равной 10, и центром в начале координат:

,    ,    .

Для метода SIA использованы следующие значения свободных параметров:

 для функции Шекеля и  для функции Растригина;

; ; ; .

Аналогично для метода HIA:

; ; ; ; ; ; .

В качестве значений свободных параметров методов CRS, MLSL, ISRES использованы их умолчательные значения, предлагаемые библиотекой NLopt [3].

Применен отличный от использованных ранее критерий окончания итераций ‑ достижение заданного числа вычислений целевой функции, равного

.

Результаты исследования представлены в таблицах 8, 9.

Таблица 8 ‑ Сравнительная эффективность методов SIA, HIA, CRS, MLSL, ISRES: функции Растригина

Таблица 9 ‑ Сравнительная эффективность методов SIA, HIA, CRS, MLSL, ISRES: функции Шекеля

Из таблиц 8, 9 следует, что предложенный метод глобальной оптимизации SIA показал очень хорошие результаты на тестовой функции Шекеля, найдя наилучшее решение среди всех рассмотренных методов. Для функции Растригина по всем критериям метод показал третий результат.

Сильное влияние на эффективность метода оказывает свободный параметр . Для функции Растригина оптимальное значение этого параметра равно 80, а для функции Шекеля ‑ 10.

Заключение

В работе предложен оригинальный гибридный метод глобальной оптимизации SIA, основанный на технологии искусственных иммунных систем. Метод представляет собой глубокую модификацию известного гибридного алгоритма HIA, который для локального поиска использует процедуру мутации. В методе SIA для этой цели применен алгоритм SUBPLEX, реализованный в известной программной библиотеке NLopt.

Разработано программное обеспечение, реализующее метод SIA, которое представляет собой многофункциональную расширяемую систему для тестирования и исследования эффективности методов оптимизации.

На тестовых функциях Розенброка, Растригина, и Химмельблау, Швефеля и Шекеля выполнено исследование эффективности алгоритма SUBPLEX, метода SIA, а также сравнительное исследование эффективности методов HIA, SIA и алгоритмов глобальной оптимизации CRS, MLSL, ISRES библиотеки NLopt. Результаты исследования показывают высокую эффективность метода SIA для целевых функций, имеющих размерность, не превышающую 50, что является достаточным для широкого круга практически важных задач оптимизации.

Для повышения эффективности метода SIA при поиске минимума функции, характер ландшафта которой неизвестен (что является типичной для приложений ситуацией), предполагается модификация метода на основе самоадаптации его ключевых свободных параметров [12]. Предполагается также разработка параллельных вариантов метода для кластерных вычислительных систем, а также для графических процессорных устройств.

Литература

1. Lucinska M., Wierzchon S.T. Hybrid Immune Algorithm for Multimodal Function Optimization // Recent Advances in Intelligent Information Systems, 2009, pp. 301-313. (http://iis.ipipan.waw.pl/2009/proceedings/iis09-30.pdf).

2. Rowan T.H. Functional Stability Analysis of Numerical Algorithms, Ph.D. Thesis, Department of Computer Sciences, University of Texas at Austin, 1990. (http://reference.kfupm.edu.sa/content/f/u/functional_stability_analysis_of_numeric_1308737.pdf).

3. The NLopt nonlinear-optimization package.

(http://ab-initio.mit.edu/wiki/index.php/NLopt).

4. Molga M., Smutnicki C. Test functions for optimization needs. 2005. (http://www.zsd.ict.pwr.wroc.pl/files/docs/functions.pdf).

5. De Castro L.N., Von Zuben F.J. Learning and Optimization Using the Clonal Selection Principle // IEEE Transactions on Evolutionary Computation, 2002, vol. 6, No. 3, pp. 231-259.

(http://www.idi.ntnu.no/emner/tdt4/2010%20articles/ clonalselection-13.10.pdf).

6. De Castro L.N., Timmis J. An Artificial Immune Network for Multimodal Function Optimization // Proceedings of IEEE Congress on Evolutionary Computation (CEC'02), 2002, May, Hawaii, vol. 1, pp. 699-674.

(http://neuro.bstu.by/ai/To-dom/My_research/Papers-0/Etc-done/wais02.pdf).

7. Yildiz A.R. A novel hybrid immune algorithm for global optimization in design and manufacturing // Robotics and Computer-Integrated Manufacturing, 2009, vol. 25, pp. 261-270.

8. Kelsey J., Timmis J. Immune Inspired Somatic Contiguous Hypermutation for Function Optimisation, 2003.

(http://www.cs.york.ac.uk/rts/docs/ GECCO_2003/papers/2723/ 27230207.pdf).

9. Agiza H.N., Hassan A.E., Salah A.M. An Improved Version of opt-AiNet Algorithm (I-opt-AiNet) for Function Optimization // IJCSNS International Journal of Computer Science and Network Security, March 2011, vol. 11, No. 3, pp. 80-85. (http://paper.ijcsns.org/07_book/201103/20110313.pdf).

10. Aragon V.S., Escuivel S.C., Coello Coello C.A. Solving Constrained Optimization using a T-Cell Artificial Immune System, 2008.

(http://polar.lsi.uned.es/revista/index.php/ia/article/viewFile/579/563).

11. Nelder J. A., Mead R. A simplex method for function minimization // The Computer Journal, 1965, Vol. 7, p. 308-313.

12. Eiben A. E., Michalewicz Z., Schoenauer M., Smith J. E. Parameter Control in Evolutionary Algorithms / Parameter Setting in Evolutionary Algorithms.- Springer Verlag, 2007, pp. 19-46.