Оптимізація роєм частинок (PSO) — Роїна інтелект

1. Походження та інтуїція

Ідея Кеннеді та Еберхарта надзвичайно проста: хороші рішення притягують сусідів. Так само як птах коригує свій політ, аби тримати темп зі зграєю, і водночас пам'ятає місця, де сам знаходив найкращу їжу, кожна частинка PSO несе соціальну пам'ять (найкраща відома позиція рою — gbest) і когнітивну пам'ять (найкраща власна позиція — pbest).

Біологічна аналогія йде ще далі. У реальній зграї жоден птах не має глобального знання — кожен бачить лише найближчих сусідів. Стандартний PSO абстрагує це у повністю зв'язану соціальну топологію, де кожна частинка може бачити gbest, але локальні варіанти топології обмежують зв'язок найближчими сусідами, породжуючи багатшу динаміку дослідження.

Історична довідка: Кеннеді та Еберхарт представили PSO на ICNN-95 (IEEE International Conference on Neural Networks). Алгоритм надихався одночасно моделями зграй (Boids Рейнольдса) і теорією соціального порівняння психолога Леона Фестінгера — ідеєю, що індивіди коригують переконання, порівнюючи себе з однолітками.

PSO належить до широкої сім'ї алгоритмів роїного інтелекту поряд з Мурашиною колонією (ACO), Штучною бджолиною колонією (ABC) та Алгоритмом світляків. Усі вони мають спільну ключову властивість: близькооптимальна глобальна поведінка виникає з простих локальних правил без централізованого управління.

2. Алгоритм PSO

PSO обшукує n-вимірний неперервний простір. Рій з N частинок ініціалізується випадковим чином у межах пошuku. На кожній ітерації кожна частинка оновлює свою швидкість і позицію:

Ініціалізація: N частинок з випадковими позиціями xᵢ ∈ [x_min, x_max]ⁿ і швидкостями vᵢ ∈ [−v_max, v_max]ⁿ.
Обчислення цільової функції f(xᵢ) для кожної частинки.
Оновлення персонального найкращого: якщо f(xᵢ) < f(pbest_i), встановити pbest_i = xᵢ.
Оновлення глобального найкращого: якщо будь-який pbest_i кращий за gbest, встановити gbest = pbest_i.
Оновлення швидкості за рівнянням оновлення швидкості (Розділ 3).
Оновлення позиції: xᵢ = xᵢ + vᵢ.
Обмеження: якщо xᵢ виходить за [x_min, x_max], застосувати обробку межі.
Повторювати кроки 2–7 до збіжності або досягнення максимальної кількості ітерацій.

Алгоритм передбачає мінімізацію за угодою (для максимізації — інвертувати цільову функцію). Межі пошуку та ліміт швидкості v_max — основні специфічні для задачі налаштування.

Стан частинки

Кожна частинка i несе чотири вектори:

xᵢ — поточна позиція в n-вимірному просторі
vᵢ — поточна швидкість (напрям і швидкість руху)
pbest_i — найкраща позиція, знайдена частинкою i за весь час
f(pbest_i) — значення цільової функції у pbest_i

Рій також підтримує єдину спільну позицію gbest — найкращу позицію, знайдену будь-якою частинкою за весь час.

3. Рівняння оновлення швидкості

Ядром PSO є правило оновлення швидкості. Канонічна форма з інерційною вагою (Shi та Eberhart, 1998):

vᵢ(t+1) = ω \cdot vᵢ(t) + c₁ \cdot r₁ \cdot (pbest i - xᵢ(t)) + c₂ \cdot r₂ \cdot (gbest - xᵢ(t))

де r₁ і r₂ — незалежні рівномірно розподілені випадкові числа з [0, 1], що перевибираються на кожному кроці для кожного виміру. Три члени відповідають трьом окремим силам:

ω·v

Інерція

Імпульс: частинка продовжує рухатись у поточному напрямку. ω > 1 сприяє дослідженню; ω < 1 — збіжності.

c₁·r₁

Когнітивна складова

Частинка притягується до власної найкращої позиції — "ностальгія" за минулим особистим успіхом.

c₂·r₂

Соціальна складова

Частинка притягується до глобального найкращого рою — "слідування натовпу" до найкращої відомої зони.

Після обчислення нової швидкості позиція оновлюється просто:

xᵢ(t+1) = xᵢ(t) + vᵢ(t+1)

Обмеження швидкості

Без обмежень на швидкість частинки можуть летіти довільно швидко і перестрибувати оптимум. Обмеження швидкості стримує кожен вимір:

vᵢ,d = clamp(vᵢ,d, -v max,d, +v max,d)

Поширена евристика: v_max,d = (x_max,d − x_min,d) / 2 — половина діапазону пошуку по кожному виміру. Це запобігає передчасній збіжності через надмірну швидкість дослідження.

Стохастичні компоненти

Випадкові масштаби r₁ і r₂ перевибираються для кожної частинки, кожного виміру, на кожній ітерації. Це суттєво: без випадковості частинки у симетричній конфігурації могли б заблокуватися в ідентичних патернах руху і рій перетворився б на детерміноване оновлення у стилі градієнта. Стохастичність також забезпечує неявні перезапуски, коли частинка далека як від pbest, так і від gbest.

4. Ключові параметри

Параметр	Символ	Типове значення	Вплив
Розмір рою	N	20–50	Більший N → краще покриття, більше обчислень на ітерацію
Інерційна вага	ω	0.4–0.9 (лінійний спад)	Висока ω досліджує; низька ω концентрується. Спад від 0.9→0.4 — найпоширеніший графік.
Когнітивний коефіцієнт	c₁	1.5–2.0	Вага притягання до персонального найкращого. Занадто висока → частинки розходяться незалежно.
Соціальний коефіцієнт	c₂	1.5–2.0	Вага притягання до глобального найкращого. Занадто висока → передчасна збіжність до локальних оптимумів.
Ліміт швидкості	v_max	½·(x_max−x_min)	Запобігає "вибуху" частинок. Занадто малий → повільна збіжність.
Макс. ітерацій	T	200–1000	Залежить від задачі; відстежуйте покращення gbest на ітерацію.

Графіки інерційної ваги

Інерційна вага ω — найбільш чутливий параметр. Фіксована низька ω (≈ 0.4) спричиняє швидку, але можливо передчасну збіжність. Лінійно спадна інерційна вага (LDIW) балансує дослідження та використання:

ω(t) = ω max - (ω max - ω min) \cdot (t / T) // Типово: ω_max = 0.9, ω_min = 0.4, T = максимальна кількість ітерацій

Нелінійні графіки (експоненційний спад, хаотичні відображення, адаптивна ω на основі різноманітності рою) можуть перевершити LDIW на специфічних задачах, але додають складності.

5. Збіжність та стабільність

Clerc та Kennedy (2002) формалізували умови збіжності для оригінального PSO. Ігноруючи обмеження швидкості та розглядаючи оновлення як стохастичне різницеве рівняння, траєкторія окремої частинки в очікуванні збігається тоді і тільки тоді, коли:

// Необхідна і достатня умова збіжності: 0 < ω < 1 0 < c₁ + c₂ < 4(1 + ω) / (1 - ω²) // спрощена межа // На практиці: ω = 0.729, c₁ + c₂ \leq 4 (значення скорочення Clerc)

Без контролю параметрів PSO осцилює і може не збігатись. Аналіз показує, що оптимальна область у просторі параметрів (ω, c₁, c₂) утворює трикутну область, відому як трикутник збіжності.

Стагнація та передчасна збіжність

Найпоширеніший режим збою — передчасна збіжність: весь рій збігається до локального оптимуму до дослідження повного простору пошуку. Це відбувається, коли:

Соціальний коефіцієнт c₂ занадто великий відносно c₁, змушуючи частинки поспішати до поточного gbest.
Рою бракує різноманітності — всі pbest кластеризуються поблизу одного регіону.
Ландшафт цільової функції має оманливий локальний оптимум, який притягує рій.

Стратегії пом'якшення включають випадковий перезапуск застиглих частинок, локальну топологію (lbest PSO, Розділ 6) або збільшення розміру рою N.

6. Варіанти: CPSO, APSO, lbest

PSO з коефіцієнтом скорочення (CPSO)

Clerc (1999) вивів коефіцієнт скорочення χ, який гарантує збіжність без явного обмеження швидкості:

φ = c₁ + c₂, φ > 4 χ = 2 / |2 - φ - sqrt(φ² - 4φ)| // Стандартні значення: c₁ = c₂ = 2.05 \to φ = 4.1 \to χ \approx 0.7298 vᵢ(t+1) = χ \cdot [vᵢ(t) + c₁\cdotr₁\cdot(pbest-xᵢ) + c₂\cdotr₂\cdot(gbest-xᵢ)]

CPSO з χ ≈ 0.73 та c₁ = c₂ = 2.05 математично еквівалентний PSO з інерційною вагою ω = 0.729 та c₁ = c₂ = 1.495. При однакових випадкових насінинах обидві форми дають ідентичні траєкторії.

PSO з локальним найкращим (lbest)

Замість єдиного спільного gbest, кожна частинка i спілкується лише з k найближчими сусідами (зазвичай k = 2 або 3 у кільцевій топології). Кожна частинка використовує lbest — найкращу позицію, знайдену у своєму локальному сусідстві:

vᵢ(t+1) = ω\cdotvᵢ + c₁\cdotr₁\cdot(pbest i -xᵢ) + c₂\cdotr₂\cdot(lbest i -xᵢ)

lbest PSO збігається повільніше, але суттєво стійкіший до передчасної збіжності на мультимодальних функціях. Хороші розв'язки поширюються кільцем поступово, дозволяючи різним частинам рою досліджувати різні зони притягання одночасно.

Адаптивний PSO (APSO)

APSO (Zhan та ін., 2009) відстежує еволюційний стан рою (дослідження / використання / збіжність / вискакування) за допомогою міри рангу пристосованості та щільності частинок і адаптує ω та c₁, c₂ у реальному часі. У стані дослідження ω збільшується; у стані збіжності ω зменшується і соціальна вага c₂ підвищується. Це усуває потребу в офлайн-підборі параметрів.

PSO на кістках (Bare-Bones PSO)

Kennedy (2003) запропонував версію з гаусівською вибіркою, що повністю усуває швидкість:

μᵢ = (pbest i + gbest) / 2 σᵢ = |pbest i - gbest| xᵢ(t+1) = Normal(μᵢ, σᵢ) // нова позиція вибирається безпосередньо

Bare-Bones PSO дивовижно конкурентоспроможний на багатьох тестових задачах і не потребує жодного підбору параметрів.

7. PSO проти GA та SA

Властивість	PSO	Генетичний алгоритм	Імітація відпалу
Популяція	Так (N частинок)	Так (N особин)	Ні (один розв'язок)
Потрібен градієнт	Ні	Ні	Ні
Використовувана пам'ять	pbest, gbest	Немає (тиск відбору)	Немає
Механізм дослідження	Інерція + випадкові компоненти	Схрещування + мутація	Випадкове збурення, графік температури
Дискретні простори	Можливо (BPSO)	Природньо (бінарний GA)	Природньо (SA на комбінаторних)
Складність реалізації	Дуже низька	Середня	Низька
Швидкість збіжності	Швидка (десятки ітерацій)	Середня	Повільна (багато збурень)
Гарантія глобального оптимуму	Ні	Ні	Так (при T → 0 достатньо повільно)

На практиці PSO часто надається перевага для неперервної оптимізації малої та середньої розмірності (≤ 100 вимірів) завдяки простій реалізації, швидкій збіжності та малій кількості гіперпараметрів. GA відмінно підходять для комбінаторних задач і задач, що вимагають структурних варіацій (наприклад, еволюція топології нейронних мереж через NEAT). SA відмінний, коли є чітка аналогія відпалу або коли потрібен один якісний розв'язок після довгого виконання.

8. Реалізація на JavaScript

// Оптимізація роєм частинок — чистий 2D приклад
// Мінімізує функцію Растригіна f(x,y) = 20 + x²−10cos(2πx) + y²−10cos(2πy)
// Глобальний мінімум в (0,0) = 0

function rastrigin(x, y) {
  return 20 + x*x - 10*Math.cos(2*Math.PI*x)
            + y*y - 10*Math.cos(2*Math.PI*y);
}

// Конфігурація PSO
const N       = 40;      // розмір рою
const DIMS    = 2;       // розмірність
const LB      = -5.12;   // нижня межа по кожному виміру
const UB      =  5.12;   // верхня межа по кожному виміру
const W_MAX   = 0.9;     // початкова інерція
const W_MIN   = 0.4;     // кінцева інерція
const C1      = 1.5;     // когнітивний коефіцієнт
const C2      = 1.5;     // соціальний коефіцієнт
const VMAX    = (UB - LB) / 2;  // ліміт швидкості
const MAX_IT  = 300;

// Ініціалізація
const pos   = Array.from({length: N}, () => Array.from({length: DIMS}, () => LB + Math.random()*(UB-LB)));
const vel   = Array.from({length: N}, () => Array.from({length: DIMS}, () => (Math.random()-0.5)*2*VMAX));
const pbest = pos.map(p => [...p]);
const pfval = pbest.map(p => rastrigin(p[0], p[1]));

let gbest = pbest[pfval.indexOf(Math.min(...pfval))].slice();
let gfval = Math.min(...pfval);

// Основний цикл
for (let t = 0; t < MAX_IT; t++) {
  const w = W_MAX - (W_MAX - W_MIN) * (t / MAX_IT);   // лінійно спадна ω

  for (let i = 0; i < N; i++) {
    for (let d = 0; d < DIMS; d++) {
      const r1 = Math.random(), r2 = Math.random();

      // Оновлення швидкості
      vel[i][d] = w * vel[i][d]
                + C1 * r1 * (pbest[i][d] - pos[i][d])
                + C2 * r2 * (gbest[d]    - pos[i][d]);

      // Обмеження швидкості
      vel[i][d] = Math.max(-VMAX, Math.min(VMAX, vel[i][d]));

      // Оновлення позиції
      pos[i][d] += vel[i][d];

      // Обробка межі: відбиття
      if (pos[i][d] < LB) { pos[i][d] = LB; vel[i][d] *= -0.5; }
      if (pos[i][d] > UB) { pos[i][d] = UB; vel[i][d] *= -0.5; }
    }

    // Обчислення
    const fval = rastrigin(pos[i][0], pos[i][1]);

    // Оновлення персонального найкращого
    if (fval < pfval[i]) {
      pfval[i] = fval;
      pbest[i] = pos[i].slice();
    }

    // Оновлення глобального найкращого
    if (fval < gfval) {
      gfval = fval;
      gbest = pos[i].slice();
    }
  }
}

console.log('Найкраща позиція:', gbest);
console.log('Найкраще значення:', gfval.toFixed(6));
// Типовий результат: значення ≈ 0.000001, позиція ≈ [0.0001, 0.0002]

Стратегії обробки меж

Чотири поширені підходи, коли частинка виходить за межі пошуку:

Поглинання: затисніть позицію до межі, зануліть компонент швидкості, що досяг межі.
Відбиття: затисніть позицію до межі, інвертуйте (і опційно масштабуйте) швидкість.
Демпфування: затисніть позицію, помножте швидкість на −λ (0 < λ < 1).
Тороїдальне: перемотайте позицію на протилежну межу — корисно для переодичних функцій.

Відбиваюча обробка меж (використана у прикладі вище) зазвичай надійна для більшості ландшафтів.

Порада з налаштування: При першому застосуванні PSO до нової задачі починайте з N = 30, ω = 0.729, c₁ = c₂ = 1.494 (значення CPSO за замовчуванням), 200 ітерацій. Якщо алгоритм рано стагнує, дещо зменшіть c₂ (до 1.3) або перейдіть до кільцевої топології (lbest). Якщо збіжність повільна, збільшіть N або додайте періодичний випадковий перезапуск для найгірших 10% частинок.

9. Застосування

Інженерне проектування

PSO широко застосовується для інженерної оптимізації, де цільові функції обчислювально дорогі, а градієнти недоступні або ненадійні. Типові приклади: оптимальні коефіцієнти ПІД-регулятора, проектування конструкцій ферм та рам (мінімізація ваги при обмеженнях на напругу), геометрія антенних решіток для мінімального рівня бічних пелюсток.

Машинне навчання та нейронні мережі

PSO може оптимізувати ваги нейронних мереж як альтернативу градієнтному зворотному поширенню, особливо коли ландшафт втрат має багато сідлових точок або коли інформація про градієнт недоступна (наприклад, недиференційовні функції активації). Мережі, навчені PSO, менш схильні потрапляти у погані локальні мінімуми на малих наборах даних.

Оптимізація гіперпараметрів

PSO розглядає гіперпараметри (швидкість навчання, регуляризацію, параметри архітектури) як виміри у неперервному просторі пошуку і оптимізує їх проти втрати крос-валідації. Природно обробляє змішані неперервні/цілочисельні виміри при поєднанні з округленням.

Багатоцільовий PSO (MOPSO)

Для задач із кількома конфліктуючими критеріями (наприклад, мінімізувати вартість і максимізувати продуктивність), MOPSO підтримує архів Парето недомінованих розв'язків. Соціальний атрактор gbest вибирається випадково з архіву зі зважуванням за щільністю архіву. З часом рій окреслює наближення до фронту Парето.

Енергетичні системи та відновлювана енергетика

Економічна диспетчеризація, оптимальний потік потужності, оптимізація розташування вітрових турбін (максимізація виробництва ферми при мінімізації втрат у кільватері) — усе це включає невипуклі, недиференційовні цільові функції, де безградієнтна природа PSO є вирішальною перевагою.

Межа масштабованості: Стандартний PSO погіршується для задач дуже великої розмірності (d > 500). Випадкова вибірка r₁ та r₂ по вимірах означає, що рій досліджує нескінченно малу частку простору. Кооперативний PSO (CPSO-H) вирішує це, розбиваючи виміри на підрої, що оптимізують свою частку, утримуючи інші при поточному значенні gbest.