Discrete optimization algorithm based on probability distribution with transformation of target values

Full Text

1. ВВЕДЕНИЕ

Решение задач оптимизации является необходимостью практически во всех сферах жизнедеятельности человека. Настоящая работа сосредоточена на комбинаторной оптимизации, которая занимается проблемами нахождения оптимума с дискретными значениями возможных решений. Одним из частных случаев этой проблемы является бинарная оптимизация, в которой элементы вектора решения могут принимать только два значения. Практическое применение таких задач весьма обширно. В области медицины бинарная оптимизация применялась для диагностики опухолей головного мозга [1], нахождения подмножеств согласованных признаков при прогнозировании эффективности реабилитации пациентов после коронавирусной инфекции [2], классификации сложных заболеваний [3], классификации аритмии по электрокардиограмме [4]. В экономической сфере для выбора издателей журналов при размещении рекламы [5], планировании рабочего процесса [6], планировании выпуска новой версии программного обеспечения [7], проектировании производственных ячеек [8]. В науке и технике бинарная оптимизация использовалась для нахождения подмножества информативных признаков при построении прогностических систем [9–11], восстановлении нагрузки в первичных распределительных сетях [12], диагностики неисправности энергосистем [13], решении проблемы позиционирования антенны [14], проектировании сварных балок [15], разделение аппаратного и программного обеспечения во встроенных системах [16]. Так же в [17] отмечается, что помимо чисто комбинаторных задач, задачи с вещественными числами могут быть представлены в двоичном виде и решены в дискретном числовом пространстве.

Для решения задач бинарной оптимизации применяют два типа методов. Первый тип – это традиционные детерминированные методы оптимизации, а второй тип основан на стохастических, недетерминированных алгоритмах. Традиционными являются методы релаксации, Лагранжа, ветвей и границ, целочисленное программирование [18, 19]. Эти методы являются трудозатратными и предназначены для решения задач небольших размерностей, что на практике встречается очень редко. Кроме того, большинство традиционных методов требуют аналитическое задание целевой функции, а также ее дифференцируемость и непрерывность. Задачи большой размерности с множеством локальных оптимумов значительно ухудшают поиск традиционных методов.

Недетерминированные методы, представляемые метаэвристическими алгоритмами [20–23], устраняют вышеперечисленные проблемы. В отличие от традиционных методов данные алгоритмы не подвергнуты “застреванию” в локальных оптимумах, в меньшей степени зависят от исходных отправных точек, не ограничены видом целевой функции и способны решать оптимизационную проблему “черного ящика” [24].

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

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

Основной научный вклад работы представлен следующими пунктами.

1. Разработан новый популяционный метаэвристический алгоритм оптимизации для поиска в дискретном пространстве. Алгоритм использует распределения вероятностей для выбора значений переменных. Распределения формируются с помощью трансформации целевых значений решений в весовые коэффициенты.
2. Эмпирически доказана эффективность предложенного алгоритма для поиска в бинарном пространстве с помощью критериев сходимости и стабильности. Статистический тест Уилкоксона показал значимое преимущество предлагаемого алгоритма по сравнению с генетическим алгоритмом и алгоритмом роящихся частиц для оптимизации унимодальных и мультимодальных тестовых функций.

Остальная часть статьи организована следующим образом. В п. 2 рассмотрены подходы и методы решения задач бинарной оптимизации с помощью метаэвристических алгоритмов; в п. 3 представлен новый алгоритм и детали его работы; в п. 4 описана экспериментальная часть исследования; в п. 5 обсуждены полученные результаты; в заключении сделаны выводы о проделанной работе.

2. БЛИЗКИЕ РАБОТЫ ПО ТЕМЕ ИССЛЕДОВАНИЯ

Наиболее популярные алгоритмы бинарной оптимизации относятся к алгоритмам роевого интеллекта. Подобно эволюционным они основаны на механизмах природы и представляют собой модель скоординированного поведения объектов, которые могут быть представителями флоры, фауны или физическими объектами. Эволюционные вычисления основаны на конкуренции и естественном отборе, тогда как роевой интеллект опирается главным образом на сотрудничество агентов [26].

Большинство алгоритмов роевого интеллекта разработаны для непрерывной оптимизации и для того чтобы осуществлять поиск в бинарном пространстве применяются механизмы адаптации, называемые бинаризацией [27]. Самым популярным методом бинаризации является использование трансформационных функций, которые переводят непрерывные значения элементов векторов решений в значения из диапазона [0, 1]. Затем применяется правило бинаризации, при котором решение преобразуется в бинарное значение из множества {0, 1}. С помощью функций трансформации были адаптированы алгоритмы роящихся частиц [17, 28], искусственных водорослей [29], шимпанзе [30], роя сальпов [31], стаи китов [32]. В [28] были исследованы различные варианты функций трансформации для алгоритма роящихся частиц. Лучшая сходимость была достигнута алгоритмом с V-образной функцией трансформации.

Метод бинаризации на основе модификации алгебраических операций преобразует вещественные операторы, используемых в формулах перемещения частиц, в их логические аналоги, что позволяет оперировать бинарными решениями. Например, вместо сложения используется операция дизъюнкции, а вместо умножения – конъюнкция. С помощью данного метода были адаптированы алгоритмы роящихся частиц [33], мозгового штурма [13], роста деревьев [34], летучих мышей [33], непрерывной муравьиной колонии [5], кукушкин поиск [35], черной дыры [36].

Квантовый метод бинаризации тоже преобразует операторы непрерывного алгоритма. В этом методе каждое допустимое решение имеет позицию и вектор квантования, который содержит вероятности принять значение 1 для соответствующего элемента решения. Вектор квантования обновляются с учетом положений глобальных и локальных лидеров. Используя данный метод, были адаптированы алгоритмы роящихся частиц [37], искусственных водорослей [7], гравитационный поиск [38], муравьиной колонии [39].

Среди алгоритмов эволюционного интеллекта для бинарной оптимизации широко применялся генетический алгоритм [40–42]. Также были предложены его модификации, так, например, в работе [43] представлен гибрид на основе генетического алгоритма и алгоритма роящихся частиц. Сначала оба алгоритма находят решения независимо друг от друга, а затем результаты объединяются с помощью метода средневзвешенной комбинации. После этого применяется локальный поиск для нахождения окончательного решения.

Оценка эффективности алгоритмов в большинстве исследований проводилась при решении определенных прикладных задач. Для объективной оценки работы алгоритмов применяют тестовые функции, которые позволяют определить эффективность при нахождении оптимума различных целевых функций, например, унимодальных, мультимодальных, овражных, разрывных, выпуклых, вогнутых. При бинарной оптимизации применяют тестовые функции для поиска в непрерывном пространстве. Бинарное пространство поиска образуют путем дискретизации непрерывного и последующим бинарным кодировании дискретных значений [28, 44–46].

3. НОВЫЙ ДИСКРЕТНЫЙ АЛГОРИТМ ОПТИМИЗАЦИИ

В настоящей работе рассматривается проблема оптимизации, в которой минимизируется критерий эффективности. В данном разделе представлен оригинальный дискретный метаэвристический алгоритм оптимизации на основе распределения вероятностей с трансформацией целевых значений (Probability Distributions with Transformation of target values, PDT). Алгоритм является итерационным, где на каждой итерации формируется вероятностная модель. Модель определяет вероятность появления конкретного дискретного значения переменной. Вероятности формируются на основе частоты появления дискретного значения каждой переменной среди решений популяции, причем меньшее значение целевой функции должно увеличивать вклад решения в повышение вероятности. Для этого вводятся трансформационные функции, которые переводят значение целевой функции решения популяции в весовой коэффициент. На рис. 1 представлена блок-схема алгоритма.

Рис. 1. Блок-схема алгоритма.

На этапе инициализации определяется начальная популяция решений Pop случайным или иным образом. Далее рассчитываются весовые коэффициенты решений популяции w. Весовой коэффициент принимает значение из диапазона [0, 1], чем меньше целевое значение решения, тем больше значение w. Для того чтобы сформировать весовые коэффициенты из целевых значений предлагается использовать функции трансформации. В качестве таких функций, например, могут быть использованы следующие: T_L – линейная функция, T_S – сигмоида, T_Q – квадратичная функция, T_Th – гиперболический тангенс. Графики функций представлены на рис. 2. Область определения функций ограничена отрезком [f_min, f_max], где f_min и f_max – минимальное и максимальное значение целевой функции в популяции соответственно, а f – текущее значение. Аналитические выражения функций трансформации показаны ниже:

$T_L\left(f\right)=\frac{f_{\max }-f}{f_{\max }-f_{\min }}$,

$T_s\left(f\right)=\frac{1}{1+e^{\left(10\left(f-f_{\min }\right)/\left(f_{\max }-f_{\min }\right)-5\right)}}$,

$T_Q\left(f\right)={\left(\frac{f-f_{\min }}{f_{\max }-f_{\min }}\right)}^2$,

$T_{Th}\left(f\right)=\left|th\left(4\left(\frac{f-f_{\min }}{f_{\max }-f_{\min }}-1\right)\right)\right|$.

Рис. 2. Графики функций трансформаций для перевода целевых значений в веса решений.

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

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

Ниже представлено пошаговое описание алгоритма.

Вход: Установить размер популяции N, число итераций MaxIter, вероятность мутации p_a и функцию трансформации T. Обозначим l_jk k-е дискретное значение j-й переменной.

Выход: R – найденное решение.

Начало

Шаг 1. Инициализация.

Случайным или иным образом сгенерировать популяцию решений Pop = [Pop₁, Pop₂, …, Pop_N] и вычислить соответствующее целевое значение f = [f₁, f₂, …, f_N].

Шаг 2. Инициализировать счетчик итераций t = 1. Начало итерационного процесса.

Шаг 3. Вычислить весовые коэффициенты решений популяции Pop с помощью функции трансформации T :

$w_i=T\left(f_i\right),$

где i = 1, ..., N.

Шаг 4. Вычислить распределения вероятностей.

Для каждого дискретного значения k переменной j определить сумму весовых коэффициентов решений Pop, которые принимают данное дискретное значение k. Обозначим такую сумму S_jk, где j = 1, …, n, k = 1, …, m:

$S_{jk}=\sum _{i=1}^N{w}_i\cdot \left\{\begin{array}{l}\mathrm{1,}еслиPo{p}_{ij}=l_{jk}\\ \mathrm{0,}иначе\end{array}\right.$.

Вычислить эмпирическую вероятность появления k-го значения переменной j:

$P_{jk}=\frac{S_{jk}}{{\displaystyle \sum _{k=1}^m{S}_{jk}}}$.

Шаг 5. Сгенерировать новые решения.

Формируется популяция новых решений Pop^new на основе вероятностей P каждой переменной:

$Po{p}_{ij}^{\mathrm{new}}=l_{jk}$,

где k удовлетворяет условию p_k–1 < rand (0,1) ≤ p_k,

$p=\left[\mathrm{0,}{P}_{j1},P_{j1}+P_{j2},\dots ,{\sum }_{k=1}^m{P}_{jm}\right]$,

j = 1, …, n, k = 1, …, m.

Шаг 6. Изменить новые решения (Мутация).

Изменить значения элементов векторов решений Pop^new с вероятностью p_a. Новые значения выбираются случайным образом из области значений элемента вектора решений.

Шаг 7. Вычислить значение целевых функций f_i^new, где i = 1, …, N, для каждого решения Pop^new.

Шаг 8. Обновить популяцию Pop путем выбора N лучших решений из множества решений Pop ∪ Pop^new. Обновить значения целевых функций f_i согласно решениям Pop.

Шаг 9. Проверка остановки алгоритма.

Если t < MaxIter, то t = t + 1 и перейти на Шаг 3, иначе перейти на Шаг 10.

Шаг 10. Извлечь в R лучшее решение из популяции Pop.

Конец

4. ЭКСПЕРИМЕНТ

В настоящем разделе представлены эксперименты c разработанным дискретным алгоритмом оптимизации на основе распределения вероятностей с трансформацией целевых значений. Алгоритм тестировался для бинарной проблемы оптимизации, т. е. когда переменные принимают только два значения. В экспериментальном исследовании использовались восемнадцать различных унимодальных и мультимодальных эталонных функции, широко применяемых для тестирования алгоритмов оптимизации [28, 44–46]. В табл. 1 представлены их характеристики, а на рис. 3 графики в двумерном пространстве поиска. Функции f₁–f₁₁ являются унимодальными, т. е. содержат только один глобальный оптимум. Функции f₁₂–f₁₈ являются мультимодальными и содержат один глобальный и множество локальных оптимумов, число которых экспоненциально растет с увеличением размерности задачи. Эксперимент проводился согласно методике работы [44].

Рис. 3. Графики тестовых функций в двумерном пространстве поиска.

Таблица 1. Тестовые функции эксперимента

Реализация алгоритма осуществлялась на языке MATLAB в среде программирования MATLAB R2022b. Программа доступна по ссылке https://cloud.tusur.ru/index.php/s/395znYyx87rRoDP. Эксперимент проводился на персональном компьютере под управлением операционной системы Windows 10 с 8 Гб оперативной памяти и процессором Intel Core i7-12700.

4.1. Дискретизация непрерывных значений

Поскольку алгоритм является дискретным и оперирует в эксперименте бинарными векторами решений, элементы которых принимают значение 0 или 1, проводится кодировка вещественных значений бинарным вектором. Процедура перевода бинарного вектора решения в значения вещественных переменных представлена на рис. 4. Данная процедура выполняется всякий раз, когда алгоритму необходимо рассчитать значение целевой функции. Количество переменных в эксперименте 5, количество битов для кодирования значения каждой переменной – 15. Таким образом, величина бинарного вектора решений составляет n = 5 × 15 = 75 элементов. Количество дискретных значений, которое может иметь каждая переменная, соответствует 2¹⁵. Эти значения определяются с помощью равномерного квантования на диапазоне поиска переменной. Шаг дискретизации определяется следующим образом:

$\Delta h=\frac{R_{\max }-R_{\min }}{2^{15}-1}$,

где R_min и R_max – нижняя и верхняя граница диапазона значений переменной соответственно. Фактически, бинарное значение переменной – это бинарное представление порядкового номера дискретного значения на диапазоне [R_min, R_max] с шагом дискретизации Δh.

Рис. 4. Перевод бинарного вектора решения в непрерывный вектор для вычисления значения целевой функции.

Если переменная x кодируется бинарным вектором [b₁, …, b₁₅], то вещественное значение этой переменной определяется следующим образом:

$x=R_{\min }+\Delta h\sum _{i=1}^{15}{b}_i\cdot 2^{i-1}$.

4.2. Критерии эффективности

Для оценки эффективности работы алгоритма применялись два критерия [47]. Первый оценивает сходимость алгоритма и определяется средним отклонением найденного целевого значения от фактического:

$E=\frac{1}{n_{run}}\sum _{i=1}^{n_{run}}\left|f_i-f^{\text{'}}\right|$,

где n_run – количество запусков алгоритма; f_i – найденное алгоритмом значение целевой функции в i-м запуске; f ʹ – фактическое значения оптимума целевой функции. Второй критерий оценивает стабильность работы недетерминированного алгоритма и определяется среднеквадратичным отклонением найденного оптимума целевой функции:

$STD=\sqrt{\frac{1}{n_{run}}\sum _{i=1}^{n_{run}}{\left(f_i-M\right)}^2}$,

где M – среднее значение целевой функции;

$M=\frac{{\displaystyle \sum _{i=1}^{n_{run}}{f}_i}}{n_{run}}$.

Меньшее значение обоих критериев соответствует лучшему значению эффективности.

Кроме вышеприведенных критериев в работе представлены графики сходимости алгоритмов, позволяющие оценить скорость сходимости стохастических алгоритмов и показывающие зависимость критерия Е от итерации [28, 44, 45]. Значение E, приводимое на графиках, является средним значением по запускам алгоритма.

4.3. Выбор функции трансформации целевых значений

Для выбора функции трансформации целевых значений в весовые коэффициенты решений были использованы следующее функции: T_L – линейная функция, T_S – сигмоида, T_Q – квадратичная функция, T_S – гиперболический тангенс.

Алгоритм PDT с разными функциями трансформации был использован для поиска оптимума тестовых функций. Было осуществлено 30 запусков на каждой тестовой функции. Полученные значения критериев эффективности приведены в табл. 2.

Таблица 2. Оценка эффективности алгоритма PDT с различными функциями трансформации

Для улучшения оценки эффективности эволюционных алгоритмов в [47] отмечается, что необходимо проводить статистические тесты. Недостаточно сравнивать алгоритмы по значениям E и STD [48], необходимо провести статистический тест, чтобы доказать, что предлагаемый новый алгоритм представляет собой значительное улучшение по сравнению с другими существующими методами.

Чтобы судить о том отличаются ли статистически значимо результаты работы алгоритма с различными функциями трансформации друг от друга, был проведен непараметрический статистический тест Фридмана на уровне значимости α = 0.05. Значения асимптотической значимости p-value, которые меньше 0.05, можно рассматривать как убедительное свидетельство против нулевой гипотезы H₀ [47].

Тест Фридмана множественных сравнений не выявил отклонение гипотезы H₀ для обоих критериев. Гипотеза H₀ здесь утверждение об отсутствии значимых различий между вариантами алгоритма с различными функциями трансформации. Асимптотическая значимость для критерия Е соответствует значению p-value = 0.757, а для критерия STD значению p-value = 0.590. Таким образом, выбор рассмотренных функций трансформации существенно не повлияет на эффективность работы алгоритма. В дальнейшем будет использоваться линейная функция.

4.4. Параметры эксперимента

Эффективность предлагаемого алгоритма PDT оценивалась в сравнении с такими популярными алгоритмами оптимизации как генетический алгоритм (GA) и бинарный алгоритм роящихся частиц (BPSO). Алгоритмы выполнялись в одинаковых условиях. Общие настройки имели следующие значения. Размер популяции – 30, количество итераций – 100, количество переменных – 5, число бит на одну переменную – 15, количество запусков алгоритма на каждую тестовую функцию – 30. Специфичные параметры алгоритмов GA и BPSO были установлены в значения, рекомендованные в [28, 45]. Значения специфичных параметров приведены в табл. 3.

Таблица 3. Значения параметров алгоритмов

4.5. Результаты эксперимента

В результате выполнения эксперимента были получены значения критериев эффективности каждого алгоритма. Данные значения приведены в табл. 4. Последняя строка таблицы содержит средние значения показателей. На рис. 5 показаны графики сходимости алгоритмов, позволяющие оценить скорость сходимости. Графики представлены в логарифмической шкале по оси критерия сходимости, что позволяет более четко отследить скорость сходимости алгоритмов на протяжении всей их работы.

Таблица 4. Оценки эффективности алгоритмов

5. ОБСУЖДЕНИЕ РЕЗУЛЬТАТОВ

Чтобы определить, значительно ли отличаются результаты эффективности предлагаемого алгоритма от аналогов, воспользуемся парным статистическим тестом Уилкоксона [47]. Нулевая гипотеза Н₀ теста утверждает отсутствие значимых различий в оценках эффективности сравниваемых алгоритмов. В табл. 5 представлены результаты сравнения. Асимптотическая значимость для критериев Е и STD при сравнении с генетическим алгоритмом и алгоритмом роящихся частиц соответствует значению p-value < 0.01. Сумма отрицательных рангов теста превалирует над положительными. Это говорит о том, что значения критериев алгоритма PDT статистически значимо меньше алгоритмов GA и BPSO на уровне значимости α = 0.01.

Таблица 5. Статистическое сравнение алгоритмов

Анализ рис. 5 показывает, что на начальных итерациях скорость алгоритма роящихся частиц для тестов f₆, f₇, f₉, f₁₅ и f₁₆ оказывается выше остальных алгоритмов, но начиная, примерно, с пятнадцатой итерации она спадает. В целом же алгоритм на основе распределения вероятностей с трансформацией целевых значений опережает по скорости своих конкурентов.

Рис. 5. Графики сходимости алгоритмов.

6. ЗАКЛЮЧЕНИЕ

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

7. БЛАГОДАРНОСТИ

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

ИСТОЧНИК ФИНАНСИРОВАНИЯ

Исследование выполнено за счет гранта Российского научного фонда № 24-21-00168, https://rscf.ru/project/24-21-00168/.

About the authors

K. S. Sarin

Author for correspondence.
Email: sarin.konstantin@mail.ru
Russian Federation, Prospect Lenina 40, Tomsk, 634050

References

Supplementary files

Supplementary Files

Action

1. JATS XML

Download

2. Fig. 1. Block diagram of the algorithm.

Download (53KB)

Indexing metadata

3. Fig. 2. Graphs of transformation functions for converting target values ​​into decision weights.

Download (16KB)

Indexing metadata

4. Fig. 3. Graphs of test functions in a two-dimensional search space.

Download (86KB)

Indexing metadata

5. Fig. 4. Transformation of a binary solution vector into a continuous vector for calculating the value of the objective function.

Download (26KB)

Indexing metadata

6. Fig. 5. Graphs of convergence of algorithms.

Download (78KB)

Indexing metadata