optimization

而世之奇偉、瑰怪、非常之觀,常在於險遠,         

而人之所罕至焉,故非有志者不能至也。《王安石.游褒禪山記》

optimization

求根找到零,最佳化找到極值。找到確切的輸入數值與輸出數值,輸出數值是最大值(最小值),稱作「最佳化」。

以函數圖形表達函數的極值:最大值就是最高的地方,最小值就是最低的地方。有時候最大值、最小值會在無限遠的地方。

中學數學談過一點點多項式函數的最佳化,比如一元二次多項式函數的最佳化(求拋物線的頂點)。大學微積分課程也談過多項式函數最佳化,比如一階導數等於零。

此處則是談一般的函數的最佳化。

容易找到極值的函數類型

unimodal function:單峰函數。先遞增、再遞減。只有遞增或者只有遞減,也算是單峰函數。適合trisection method。

convex function:凸函數。隨便劃一道割線,函數曲線總是凹下去(或平貼)。適合gradient descent。

polynomial function:多項式函數。「斜率等於零」的地方是極值、鞍點。可以手工推導公式解,也可以使用最佳化演算法求解。

大家習慣改成嚴格版本、去掉等號,讓函數曲線不會平行於座標軸、不會產生無限多極值。

optimization - 前後包夾類型

概論

斜率等於零的地方(不考慮無限遠處)是極值、鞍點。

預先微分,然後求根,即是極值、鞍點的所在之處。

這類型演算法源自二元搜尋。前後包夾,逐步縮小包圍網。

這類型演算法泛稱為bracketing method。

bisection method(二分法)

預先微分,利用二分法求根,得到極值所在之處。

如果無法預先推導微分公式,那麼可以即時微分。自訂dx大小,計算df,相除得到斜率df/dx。

dx和df已經不再符合數學家的定義。重新標記為Δx和Δf。

二分法當中,只需要取得斜率正負號,不需要求得斜率。只需要求得Δf,不需要除以Δx。節省計算時間,也避免浮點數誤差。

二分法求根,使用條件是兩端點位於X軸異側,f(xL) f(xR) ≤ 0。二分法求極值,使用條件變成兩端點斜率異號f′(xL) f′(xR) ≤ 0。

二分法求最大值,使用條件是f′(xL) ≥ 0且f′(xR) ≤ 0。斜率先正再負,函數先遞增再遞減,單峰函數。

二分法求最小值,函數高低顛倒、斜率正負顛倒。

也就是說,二分法適用於f(x)是單峰函數。

二分法無法求鞍點。鞍點出現條件是f′(xL) f′(xR) ≥ 0。

trisection method(三分法)

二分法改良版。中央斜率的兩點,改成三等分的兩點。

優點:程式碼變短。不需要Δx。不需要正負號。

缺點:趨近速度變慢。二分法減少1/2,三分法只減少1/3。

UVa 1476 10385 11243 11562 11613

golden section method(黃金分割法)

三分法改良版。1:1:1改成(1-ψ):(2ψ-1):(1-ψ)。

ψ是黃金比率的倒數。ψ = 1/φ = (√̅5 - 1) / 2 = 0.618...。

為什麼是黃金比率的倒數呢?我也不知道。開心就好。

優點:趨近速度稍微改善。每次減少(1-ψ)/1。

函數推廣成多個變數

演算法推廣成4ᴺ個點,時間複雜度呈指數成長,毫無實用價值。

一種改良方式:隨便抓幾點,利用加權平均數,找到更高的點。然而目前尚未有人深入研究。等你發表論文。

ensemble averaging(總體平均)

三個臭皮匠勝過一個諸葛亮。隨機挑選多個地點,得到多種極值。如果這些極值位置們都不是正確位置,那麼這些極值位置們的加權平均數,通常更接近正確位置。權重由人類手動調整。

甚至實施多種最佳化演算法,得到多種極值。

optimization - 孤軍深入類型

概論

斜率等於零的地方(不考慮無限遠處)是極值、鞍點。

預先微分,然後求根,即是極值、鞍點的所在之處。

這類型演算法源自牛頓法、不動點遞推法。孤軍深入,逐步進逼目的地。

這類型演算法泛稱為iterative method。

牛頓法求根     xnext = x - f(x)  / f′(x)
牛頓法求極值    xnext = x - f′(x) / f″(x)
不動點遞推法求根  xnext = x + λ f(x)
不動點遞推法求極值 xnext = x + λ f′(x)

函數推廣成多個變數

由於內容較多較雜,另以專門頁面介紹。請見「multivariate function」。

以下進行簡易分類。

Newton's method(牛頓法)

預先微分,然後以牛頓法求得極值、鞍點。適用於f′(x)是凸函數。

quasi-Newton method(擬牛頓法)(類牛頓法)

牛頓法改良版。割線法。二次微分的部分,改成其他類似的東西。調整方式層出不窮,人人有一套理論:

最近走紅的演算法有BHHH、BFGS。

fixed point iteration(不動點遞推法)
gradient descent(梯度下降法)

預先微分、乘上倍率λ、加x,然後以不動點遞推法求得極值、鞍點。適用於x + λ f′(x)是平緩函數。

nonlinear conjugate gradient method (非線性共軛梯度法)
adaptive gradient descent(自適應梯度下降法)

不動點遞推法改良版。前進方向是「本次的梯度」與「上次的前進方向」的加權平均。權重有許多種選擇,人人有一套理論:

最近走紅的演算法有Fletcher–Reeves。

註:名稱源自conjugate gradient method。請見本站文件「quadratic optimization」。

optimization - 地面勘查類型

概論

Show[ Plot3D[BesselJ[0, Norm[{x, y}]], {x, -4Pi, 4Pi}, {y, -4Pi, 4Pi}, PlotRange -> {-1, 1}, Boxed -> False, Mesh -> None, Axes -> False, PlotPoints -> 50, ColorFunction -> (ColorData["CherryTones"][Rescale[#3, {-4, 4}]] &)], ParametricPlot3D[{{-5, u, BesselJ[0, Norm[{-5, u}]]}, {u, -8, BesselJ[0, Norm[{u, -8}]]}}, {u, -4Pi, 4Pi}, PlotStyle -> {Red, Red}] ]
ParametricPlot3D[{v Cos[u], v Sin[u], 5 Cos[Pi/10] BesselJ[0, v]}, {u, 0, 4 (Pi/2)}, {v, 0, 15}, Boxed -> False, Axes -> False, Mesh -> None]

選一個起點,一步一腳印,走向極值。

這類型演算法源自iterative method與metaheuristics。

hill climbing(登山法)

沿著函數圖形的表面前進。隨便往某個方向跨出一步,確定是往上,就直直走;確定是往下,就不走。最後成功登頂。

步伐太大,會越過山峰,走往低處。只要步伐持續變小,仍可攻頂。但是步伐變小太快,便會在山腰掙扎,無法登頂。取捨兩難。

選擇各種起點,登上各個山峰。由於無法預測最高峰的位置,只好嘗試多種起點、實施多次登山法,依賴運氣尋找最高峰。嘗試越多起點,越能找到最高峰,越能避免落入局部極值。

缺點:可能停在山腰。可能只找到局部極值。

simulated annealing(模擬退火法)

登山法改良版。允許往下走,以便翻山越嶺。

隨便往某個方向跨出一步,確定是往上,就走;確定是往下,以機率exp(Δf⋅t)決定走或不走,其中Δf是上升高度(往下走時是負值),t是時刻(大於等於1)。大致上來說,往下越陡就越不想往下,年紀越大就越不想往下。

原論文只找山谷、未找山峰,因而取名為「退火」。原論文沒有時刻t,而是溫度T;溫度下降,因此機率公式是exp(Δf/T) 。

缺點:可能停在山腰。可能只找到局部極值。

UVa 10228

gradient descent(梯度下降法)

登山法改良版。直接朝著梯度方向走,不必試誤。

X座標、Y座標、……分開處理,互不干涉。當前位置,求得相鄰高度差,座標大的高度減去座標小的高度,正負值可判別相對高低。當前位置座標,加上相鄰高度差,就是往上走。若相鄰高度差越多,則前後座標差越多。坡越陡、跨越遠、走越快。

Plot[BesselJ[0, Norm[{-5, y}]], {y, -4Pi, 4Pi}, PlotRange -> {-1,1}] Plot[BesselJ[0, Norm[{x, -8}]], {x, -4Pi, 4Pi}, PlotRange -> {-1,1}]

梯度是相鄰高度差再除以dx,功效差不多。因為數學家喜歡梯度,所以採用梯度。為了得到梯度,函數必須一次可微。

梯度大小是傾斜程度。梯度方向是最陡的方向,是等高線的垂直方向。朝梯度方向走會上升、得極大值,朝梯度反方向走會下降、得極小值。步伐太大,會走之字路線,無傷大雅。

原論文只找山谷、未找山峰,因而取名為「下降」。

優點:方向明確,不必隨機亂走試誤,攻頂速度更快。缺點:可能停在山腰。可能只找到局部極值、鞍點。

註:梯度下降法可以視作不動點遞推法!

adaptive gradient descent(自適應梯度下降法)

梯度下降法改良版。隨時調整前進方向、隨時調整步伐大小。調整方式層出不窮,人人有一套理論:

最近走紅的演算法有Adagrad、Adam、momentum。

優點:自動調整步伐大小,緩急自如。缺點:計算量略增。

註:自適應梯度下降法可以視作非線性共軛梯度法!

line search(直線搜尋)

朝著指定方向走,走到最高處。

偏微分等於零的地方,就是該方向的極值、鞍點。

偏微分等於零的地方,梯度與方向恰好垂直。

如果可以推導公式解,那麼就能一步到位。

如果無法推導公式解,那麼只能步步登高。

優點:有了公式解,不必煩惱步伐大小。缺點:需要指定方向。

梯度下降法:當前點的法線方向。

直線搜尋:下一點的切線方向。

steepest descent(最陡下降法)

直線搜尋改良版,指定方向是梯度方向。

朝著梯度方向走,走到最高處。

相鄰兩步的方向,恰好互相垂直。

註:最陡下降法是梯度下降法與線性搜尋兩者並用!

coordinate descent(座標下降法)

直線搜尋改良版。指定方向是座標軸方向。

朝著座標軸方向走,走到最高處。每個座標軸不斷輪流。

優點:更容易推導公式解。缺點:故意繞遠路,攻頂速度更慢。

註:座標下降法可以視作Gauss–Seidel iteration!

coordinate gradient descent(座標梯度下降法)

朝著座標軸方向走一步。每個座標軸輪流一遍。

劣於梯度下降法,毫無可取之處。然而過程宛如「Gibbs sampling」,最佳化與取樣兩個領域的橋樑!

optimization - 空中勘查類型

概論

地面勘查經常癱在山腰、卡在山丘,只好改為飛行模式。像偵察機一樣飛來飛去,不被山峰峽谷阻擋。

這類型演算法,完全憑感覺,一個比一個唬爛,連一個數學性質都不必證明。但是我們也沒有更好的方法了。面對亂七八糟的函數,只好用亂七八糟的算法。

這類型演算法泛稱為metaheuristics。這裡有實際範例

particle swarm optimization(粒子演算法)

一、散佈N個粒子。一個粒子對應一個可行解。
  position[1...N];	// x
  solution[1...N];	// f(x)
二、以個人過去最佳解、團體過去最佳解,決定粒子的速度。
  然後移動粒子,讓粒子飛往最佳解。
  velocity[i] =
   2 * rand(0 ~ 1) * (best_position[i] - position[i]) +
   2 * rand(0 ~ 1) * (best_position[best_particle_index] - position[i]);
三、更新個人最佳解、團體最佳解。
  best_solution[i] = max(best_solution[i], solution[i])
    and record best_position[i];
  best_particle_solution = max_arg(best_solution[1...N])
    and record best_particle_index;
四、重複二三,直到解夠好。

附帶一提,這與粒子的真實行為完全無關。

bee colony optimization(蜜蜂演算法)

一、散佈N個食物。一個食物對應一個可行解。
  position[1...N];	// x
  solution[1...N];	// f(x)
二、N隻蜜蜂分頭前往N個食物並返回,但是腦中記得的位置會有偏差。
  position[i] +=
   rand(-1 ~ +1) * position[rand(1 ~ N except i)];
三、每隻蜜蜂皆從N個食物中挑一個去,機率由解的好壞決定。
  probability[i] = solution[i] / ( solution[1] + ... + solution[N] )
  返回時腦中記得的位置一樣會有偏差。公式同二。
四、如果某個食物連續K個回合沒去,食物自動消滅。
  隨機散佈食物,補足食物直到N個。
五、重複二三四,直到解夠好。

附帶一提,這與蜜蜂的真實行為完全無關。

constrained optimization

天將降大任於斯人也,必先苦其心志,勞其筋骨,餓其體膚,   

空乏其身,行拂亂其所為,所以動心忍性,曾益其所不能。《孟子》

constrained optimization

「約束最佳化」。最佳化,限制輸入範圍。

輸入只能是受規範的、被指定的數值。

maximize f(x)
 x ∈ C

「約束最佳化」。最佳化、等式與不等式組,兩者合體。

輸入必須同時滿足許多道等式與不等式。

maximize f(x)   subject to g(x) = 0, h(x) ≥ 0, ...

約束最佳化當中,函數稱作「目標函數objective function」。等式與不等式組稱作「約束條件constraint」。等式與不等式可以重新整理成函數的模樣,稱作「約束函數constraint function」。

通常會存在太過完美的輸入,所以只好抑制。

施加限制,以求得更對味的答案。

Show[ Plot3D[(Cos[(x-y)/2 - 3] + Cos[y/4])/2 - 1/2, {x, -4Pi, 4Pi}, {y, -4Pi, 4Pi}, PlotRange -> {-1, 1}, Boxed -> False, Mesh -> None, Axes -> False, PlotPoints -> 50, ColorFunction -> (ColorData["CherryTones"][Rescale[#3, {-4, 4}]] &)], ParametricPlot3D[{{u,-4Pi,0}, {-4Pi,u,0}}, {u, -4Pi, 4Pi}, PlotStyle -> {Black, Black}] ]
Show[ Plot3D[BesselJ[0, Norm[{x, y}]], {x, -4Pi, 4Pi}, {y, -4Pi, 4Pi}, PlotRange -> {-1, 1}, Boxed -> False, Mesh -> None, Axes -> False, PlotPoints -> 50, ColorFunction -> (ColorData["GrayTones"][Rescale[#3, {-6, 2}]] &)], Plot3D[BesselJ[0, Norm[{x, y}]] + 0.01, {x, -4Pi, 4Pi}, {y, -4Pi, 4Pi}, PlotRange -> {-1, 1}, Boxed -> False, Mesh -> None, Axes -> False, PlotPoints -> 50, RegionFunction -> Function[{x, y, z}, Cos[(x-y)/2 - 3] + Cos[y/4] > 1], ColorFunction -> (ColorData["CherryTones"][Rescale[#3, {-4, 4}]] &)] ]
Show[ Plot3D[0, {x, -4Pi, 4Pi}, {y, -4Pi, 4Pi}, PlotRange -> {-1, 1}, Boxed -> False, Mesh -> None, Axes -> False, PlotPoints -> 50, ColorFunction -> (ColorData["GrayTones"][Rescale[#3, {-6, 2}]] &)], Plot3D[0.01, {x, -4Pi, 4Pi}, {y, -4Pi, 4Pi}, PlotRange -> {-1, 1}, Boxed -> False, Mesh -> None, Axes -> False, PlotPoints -> 50, RegionFunction -> Function[{x, y, z}, Cos[(x-y)/2 - 3] + Cos[y/4] > 1], ColorFunction -> (ColorData["CherryTones"][Rescale[#3, {-4, 4}]] &)] ]
ContourPlot[(Cos[(x-y)/2 - 3] + Cos[y/4])/2 - 1/2 == 0, {x, -4Pi, 4Pi}, {y, -4Pi, 4Pi}]

conditional gradient method(條件梯度法)
feasible direction method(可行方向法)

梯度下降法改良版。總是朝向約束條件。

p = argmin [ ∇f(x) dot (p - x) ] = argmin [ ∇f(x) dot p ]
     p∈C                            p∈C

xnext = x + step ⋅ (p - x)

從約束條件當中,找出最符合當前梯度方向、最靠近當前位置的地點(以點積大小來判斷),然後朝該地點走一步。

找到該地點,是另一個最佳化問題。兩層最佳化。

如果約束函數很單純,則推導公式解,得到該地點。如果約束函數很複雜,則以梯度下降法,找到該地點。

projected gradient method(投影梯度法)

梯度下降法改良版。每當偏離約束條件,立即歸位。

p = x + step ⋅ ∇f(x)

xnext = argmin ‖ p - xnext ‖
       xnext∈C

朝著梯度方向走一步,然後立即垂直投影到約束條件上面。換句話說,若離開約束條件,則走最短直線抵達約束條件;若未離開約束條件,則維持現狀。

垂直投影,是另一個最佳化問題。兩層最佳化。

如果約束函數很單純,則推導垂直投影公式,一步到位。如果約束函數很複雜,則以梯度下降法慢慢走入約束條件:約束函數的梯度方向,權且作為垂直投影方向。

Lagrange multiplier(拉格朗日乘數)

當約束條件是等式,約束最佳化可以化作微分方程組求解。

f(x,y) 求極值。必須滿足 g(x,y) = 0。
湊出 L(x,y,λ) = f(x,y) + λ g(x,y)。

f(x,y) 的極值,等同 L(x,y,λ) = f(x,y) + λ g(x,y) 的極值。
欲求極值:
對 x 偏微分,讓斜率是 0。
對 y 偏微分,讓斜率是 0。

不管 λ 如何變化,λ g(x,y) 都是零,L(x,y,λ) 永遠不變。
欲求永遠不變的地方:
對 λ 偏微分,讓斜率是 0。

三道偏微分方程式聯立,其解涵蓋了(不全是)所有符合約束條件的極值。
⎧ ∂/∂x L(x,y,λ) = 0
⎨ ∂/∂y L(x,y,λ) = 0
⎩ ∂/∂λ L(x,y,λ) = 0

上述微分方程組,簡寫成「梯度等於零」、「求駐點」。
∇L(x,y,λ) = 0

L(x,y,λ)稱作「拉格朗日函數Lagrangian」。λ稱作「拉格朗日乘數Lagrange multiplier」。

微分方程組的解(拉格朗日函數的駐點),涵蓋所有正確答案,卻也可能包含錯誤答案。必須驗證答案。檢查每個答案的鄰近函數值(篩去鞍點)、比較每個答案的函數值(篩去局部極值)。

∇L(x,y,λ) = 0        梯度等於零、求駐點。

max L(x,y,λ)         ✘ 注意到,求駐點,不等價於求極值。
x,y,λ
max min L(x,y,λ)     ✘ 注意到,求駐點,不等價於求鞍點。
x,y   λ 

另外還可以進一步得到非常漂亮的數學性質:極值所在之處,目標函數與約束函數,兩者梯度方向共線,而梯度大小不明。(梯度大小設定成未知倍率λ。)

三式分別偏微分。
⎧ ∂/∂x L(x,y,λ) = 0
⎨ ∂/∂y L(x,y,λ) = 0
⎩ ∂/∂λ L(x,y,λ) = 0

展開L(x,y,λ) = f(x,y) + λ g(x,y)。
偏微分是線性函數,擁有分配律。
⎧ ∂/∂x f(x,y) + λ ⋅ ∂/∂x g(x,y) = 0
⎨ ∂/∂y f(x,y) + λ ⋅ ∂/∂y g(x,y) = 0
⎩ g(x,y) = 0

第一式、第二式合起來,簡寫成梯度,並且重新整理。
∇f(x,y) + λ ∇g(x,y) = 0
∇f(x,y) = -λ ∇g(x,y)

拉格朗日函數改為L(x,y,λ) = f(x,y) - λ g(x,y)。
改變正負號,式子比較帥。
∇f(x,y) = λ ∇g(x,y)

上述微分方程組,簡寫成「梯度比例」。
⎰ ∇f(x,y) = λ ∇g(x,y)
⎱ g(x,y) = 0

當約束條件是大量等式,則通通乘上倍率再加總。

f(x,y) 求極值。必須滿足 g₁(x,y) = 0 與 g₂(x,y) = 0 與……。
定義 L(x,y,λ₁,λ₂,...) = f(x,y) + λ₁ g₁(x,y) + λ₂ g₂(x,y) + ...

Karush–Kuhn–Tucker conditions(KKT條件)

約束條件是等式:稱作拉格朗日乘數。極值位於約束條件上面。

極值所在之處,目標函數與約束函數,兩者梯度方向共線,而梯度大小不明。(梯度大小設定成未知倍率λ。)

maximize f(x)   subject to g(x) = 0

---> solve ⎰ ∇f(x) = λ ∇g(x)
           ⎱ g(x) = 0

約束條件是不等式:稱作KKT條件。極值可能位於約束條件邊界或者內部。

邊界:目標函數與約束函數,兩者梯度方向共線,而梯度大小不明。(梯度大小設定成未知倍率μ。)

內部:目標函數的極值位置,已經滿足約束條件。約束條件無關緊要。(μ設定成零,讓約束條件失效。)

計算學家將兩種情況分開處理,邊界是微分方程組求解,內部是目標函數最佳化。數學家則是將兩種情況合併成一個方程組。

maximize f(x)   subject to h(x) ≥ 0

---> solve ⎧ ∇f(x) = μ ∇h(x)
           ⎨ μ ≤ 0
           ⎩ h(x) ≥ 0

極值位於約束條件邊界時,梯度方向的判斷方式:

一、最大值位於邊界:目標函數的梯度方向,離開邊界。

二、約束函數大於等於零:約束函數的梯度方向,進入邊界。

梯度方向恰好相反,梯度大小必須異號。(μ必須小於零。)

                      | ∇f(x) = μ ∇h(x) | ∇f(x) + μ ∇h(x) = 0
----------------------|-----------------|---------------------
max f(x) st h(x) ≥ 0  | μ ≤ 0           | μ ≥ 0
max f(x) st h(x) ≤ 0  | μ ≥ 0           | μ ≤ 0
min f(x) st h(x) ≥ 0  | μ ≥ 0           | μ ≤ 0
min f(x) st h(x) ≤ 0  | μ ≤ 0           | μ ≥ 0

約束條件是大量等式與不等式:各有倍率,通通加總。

maximize f(x)   subject to g₁(x) = 0, g₂(x) = 0, h₁(x) ≥ 0, h₂(x) ≥ 0

---> solve ⎧ ∇f(x) = λ₁ ∇g₁(x) + λ₂ ∇g₂(x) + μ₁ ∇h₁(x) + μ₂ ∇h₂(x)
           ⎪ g₁(x) = 0
           ⎪ g₂(x) = 0
           ⎨ μ₁ ≤ 0
           ⎪ h₁(x) ≥ 0
           ⎪ μ₂ ≤ 0
           ⎩ h₂(x) ≥ 0

penalty method(懲罰法)

約束函數取平方、自訂倍率,摻入目標函數。

     maximize f(x)               subject to g(x) = 0
---> maximize f(x)               subject to g(x) = 0 and ρ/2 g(x)² = 0
---> maximize f(x) + ρ/2 g(x)²   subject to g(x) = 0

函數變陡,梯度下降法更快找到極值,不易停在鞍點。

大家習慣讓倍率除以二,配合二次方。微分時,係數相消。

barrier method(屏障法)

約束函數取對數、自訂倍率、視情況變號,摻入目標函數。

     minimize f(x)                 subject to h(x) ≥ 0
---> minimize f(x)                 subject to h(x) ≥ 0 and
                                   -ρ log(h(x)) ≥ ∞ and ρ ≥ 0
---> minimize f(x) - ρ log(h(x))   subject to h(x) ≥ 0 and ρ ≥ 0

邊界附近,函數趨近正無限大,梯度下降法不易出界。

alternating direction method of multipliers(交替方向乘數法)

一、懲罰法:約束函數取平方,摻入目標函數。(可用可不用)

二、拉格朗日乘數:約束最佳化化作函數求駐點。

三、座標下降法:變數求極值,走到最高處;乘數偏微分求根,朝梯度方向走一步。

當前位置在約束條件內與外,梯度方向正與負。【尚待確認】

非負乘數,利用投影梯度法,利用max(μ,0)保持非負。

maximize f(x)   subject to g(x) = 0, h(x) ≥ 0
L(x,λ,μ) = f(x) + λ g(x) + μ h(x)   where μ ≥ 0
xnext = argmax { f(x) + λ g(x) + μ h(x) }
λnext = λ + step ⋅ g(x)
μnext = max(μ + step ⋅ h(x), 0)

範例:限制條件是多道一次不等式,併成一個一次方程組。

maximize f(x)   subject to Ax ≤ b
L(x,μ) = f(x) + μᵀ (Ax - b)   where μ ≤ 0
xnext = argmax { f(x) + μᵀ (Ax - b) }
μnext = min(μ + step ⋅ (Ax - b), 0)

範例:兩個函數相加。使用了懲罰法。

minimize f(x) + g(y)   subject to Ax + By = c
L(x,y,λ) = f(x) + g(y) + λᵀ (Ax + By - c) + (ρ/2) ‖Ax + By - c‖²
xnext = argmin L(x,y,λ)
ynext = argmin L(x,y,λ)
λnext = λ + step ⋅ (Ax + By - c)

interior point method(內點法)

一、屏障法:約束函數取對數,摻入目標函數。(可用可不用)

二、拉格朗日乘數:約束最佳化化作微分方程組。

三、牛頓法:多變數函數求根。

multiobjective optimization

魚,我所欲也;熊掌,亦我所欲也。    

二者不可得兼,舍魚而取熊掌者也。《孟子》

multiobjective optimization

「多目標最佳化」。找到一個輸入,同時最佳化多個函數。

⎧ minimize f(x)
⎨ maximize g(x)
⎩ maximize h(x)

通常不存在如此完美的輸入,所以只好折衷。

各取所長,以求得更均衡的答案。

scalarization(純量化)

多個函數,極值位置通常不相等。多個函數的加減乘除,極值位置毫無規律。我們難以制定任何策略。就連最佳化專家孟子也自認只能選一個函數來最佳化。這種情況下,大家只好尋求心靈勵志書籍、感恩師父感恩上人了。

⎧ minimize f(x)   Mencius bear paw
⎨ maximize g(x)  ------------------>  maximize g(x)
⎩ maximize h(x)

然而值得一提的是,凸函數相加之後,仍是凸函數;而且新極值位置,夾在原極值位置之間,不會歪太多。得到折衷方案:最佳化凸函數的和,找到極值位置。

⎧ minimize f(x)   f,g,h are convex
⎨ minimize g(x)  ------------------>  minimize f(x) + g(x) + h(x)
⎩ minimize h(x)

為了調整極值位置,可以添上權重。為了保持是凸函數,權重不能是負數。權重可以自由設定,沒有所謂最適當的權重。

注意到,除了權重,函數本身的斜率變化也會影響極值位置。即使權重均等,極值位置也通常不在正中央。因此這個手法是任憑感覺、碰碰運氣。

minimize α f(x) + β g(x) + γ h(x)  (α ≥ 0, β ≥ 0, γ ≥ 0)

衡量輕重、化作權重,故稱之為「純量化」。

regularization(正則化)

純量化的進階應用。追加其他的最佳化目標。

minimize f(x) + α g(x) + β h(x) + ...   (α ≥ 0, β ≥ 0, ...)

f(x)是原本的目標,g(x) h(x) ...是追加的目標。

大家視問題需要,訂立適當的函數。例如x的長度平方越小越好g(x) = ‖x‖²;物理學的做功越少越好g(x) = F x。

追加凸函數,導致函數變窄變陡,梯度下降法更快找到極值、不易停在鞍點。

大家視問題需要,訂立適當的權重。權重越小,函數的影響力就越小,類似於加權平均的概念。

求極值,權重的絕對大小不重要,考慮相對大小即可。原始函數的權重定為1,節省一個權重數值。

明定準則、施予規範,故稱之為「正則化」。

regularization ≠ KKT conditions

regularization:
minimize f(x) + α g(x) + β h(x) + ...   (α ≥ 0, β ≥ 0, ...)

KKT conditions:
∇[f(x) + λ g(x) + μ h(x) + ...] = 0     (μ ≥ 0, h(x) ≥ 0, ...)
regularization: α, β, ... are known parameters
KKT conditions: λ, μ, ... are unknown arguments

regularization與KKT conditions沒有任何關聯,只不過是數學式子有點像。別忘了KKT conditions還得聯立其他式子。

operator splitting method(算子分裂法)

一個函數拆成多個函數相加,每個函數輪流走一步。

minimize f(x)  --->  minimize f₁(x) + f₂(x) + f₃(x)
                     where f₁(x) + f₂(x) + f₃(x) = f(x)
x₁ = x₀ - ∇f₁(x)
x₂ = x₁ - ∇f₂(x)
x₃ = x₂ - ∇f₃(x)
x₄ = x₃ - ∇f₁(x)
x₅ = x₄ - ∇f₂(x)
:       :

大家習慣拆成特殊的函數,以便輕鬆求得梯度。

最佳化領域當中,大家習慣省略算子二字,稱作「分裂法」。

domain decomposition method(領域分解法)

算子分裂法的進階應用。一個多變數函數拆成多個少變數函數相加,每個函數輪流走一步。

minimize f(x,y,z)  --->  minimize f₁(x,z) + f₂(y,z)

大家習慣拆成特殊的函數,以便輕鬆求得極值。

最佳化領域當中,大家習慣省略領域二字,稱作「分解法」。

disciplined optimization

子路問:「聞斯行諸?」子曰:「有父兄在,如之何其聞斯行之?」

冉有問:「聞斯行諸?」子曰:「聞斯行之。」《論語》     

disciplined optimization

「紀律最佳化」。依照函數性質,採取適當的最佳化演算法。

經典函數(例如sgn、abs、min、max、pow、sqrt、exp、log),實施經典運算(例如加、減、乘、除、微、積、複合、反函數),產生有規律的函數。判斷函數具備哪些數學性質,實施相符的最佳化演算法。

新興理論,尚在發展當中。等你發表論文。

請見專論《Disciplined Convex Programming》

online optimization

故兵無常勢,水無常形。能因敵變化而取勝者,謂之神。《孫子》

online optimization

「在線最佳化」。函數時時變化,有如波浪,隨波逐流。

Plot3D[(Sin[x/2-y] + Sin[y/4])/2, {x, -4Pi, 4Pi}, {y, -4Pi, 4Pi}, PlotRange -> {-1, 1}, Boxed -> False, Mesh -> None, Axes -> False, PlotPoints -> 50, ColorFunction -> (ColorData["CherryTones"][Rescale[#3, {-4, 4}]] &)]

這類型演算法源自ensemble averaging。

mean shift(平均數移動)

一、均勻散布粒子。
二、以上次的極值座標為中心,
  取得鄰近範圍內所有粒子。(範圍自訂,通常是規定半徑長度。)
三、範圍內所有粒子各自求函數值。
四、範圍內所有粒子的座標,求加權平均(權重是函數值),得到新的極值座標。

另一種觀點是取樣。樣本(粒子)的密度,當作是函數值的大小。可以想成是找到粒子密度最高的地方。【待補文字】

particle filter(粒子濾波器)

一、以上次的極值座標為中心,
  散布n個粒子,呈高斯分布。(範圍自訂,範圍即變異數)
  粒子不用飛來飛去了,直接散布。
二、n個粒子各自求函數值。
三、n個粒子的座標,求加權平均(權重是函數值),得到新的極值座標。
四、或者不採用加權平均。
  想像n個粒子各有磁場(高斯分布),強弱由函數值大小決定,
  形成joint distribution(高斯混和分布)。
  下一回合依此分布散布粒子。

follow the leader(追隨領袖)

前幾次的函數相加後,找到極值位置,推定為本次的極值位置、或者本次的登山起點。找極值位置時,只找幾個重點位置。

重點位置隨便找,指引我們向前進──儘管非常唬爛,但是效果出眾,於是有人認真進行數學分析:

還有其他進階版本FTRL、FTRL-proximal。

evolutionary optimization

概論

函數輸入是組合或排列,函數輸出是數值。

演算法是隨意拼湊函數輸入,試著讓函數輸出是極值。嘗試多種拼湊方式,逐步改善答案。拼湊方式簡直是說故事比賽,故事材料來自自然規律、來自生物演化。

這類型演算法泛稱為evolutionary algorithm,隸屬於metaheuristics。

genetic algorithm(基因演算法)

靈感來自於染色體減數分裂的過程,優良的基因不斷遺傳下去,逐代演化出更適應環境的基因。基因演算法把答案比擬成染色體,把好的答案不斷分裂再結合,成為更好的答案。

附帶一提,這與基因的真實行為完全無關。

1.
[初始化]
一開始先隨便弄出幾個x。本例是四個。

	1010101010
	1011001011
	1110101011
	0010101000

2.
[fitness function]
根據問題特性,定義好壞程度。

	f(1010101010) = 678

3.
[selection]
隨便找個位置切一刀,每個x都被分成兩段。

	1010101  010
	1011001  011
	1110101  011
	0010101  000

4.
[crossover]
隨便找兩組你覺得夠優良的x,交叉相接變成新答案。
重複一直做,直到x數目跟原先一樣多。本例是四個。

	1010101 \/ 010  ->  1010101 -- 011
	1011001 /\ 011      1011001 -- 010 


	1010101011
	1011001010
	1110101010
	1010101000

5.
[mutation]
每個x都隨便找一個地方把數字改掉,也可以不改。

	1010111011
	1011001000
	1110101010
	1010101001

6.
重複3. 4. 5.,直到裡面有一個x是你滿意的,令f(x)最大的那個x。
1. 隨機產生N個x。
2. 計算fitness function。
3. 重複以下步驟,直到有一個x讓人滿意。
 甲、selection。
 乙、crossover。
 丙、mutation。
 丁、計算fitness function。

演算法的過程,大致就是如此。細部的實作方式、參數的調校方式,隨人話虎爛。

一開始的x的足夠豐富,多演化幾次,就可以得到不錯的結果。一開始的x足夠豐富,可以避免進入局部極值。mutation用於增加x的豐富性,以跳脫局部極值。

範例:0/1 knapsack problem

N個物品。選出其中幾個。

[初始化]
答案設計成N個二進位數字,
0代表對應的物品不在背包內,
1代表對應的物品在背包內。

[fitness function]
是背包內物品總價值,數值越大越好。

當「特定組合」具有深邃的影響力,造成最佳解、次佳解們都包含著同一(幾)套「特定組合」,此時就適合使用基因演算法。以0/1背包問題為例,一套好的物品組合方式可以有效節省背包空間,只要依賴幾套好的物品組合方式,就留有足夠的背包空間,能夠隨心所欲的放入其他物品;所有接近最佳解的答案,答案的其中一小段,都會是那幾套固定的物品組合方式──這種情況就非常適合使用基因演算法。

UVa 10715

範例:minimum spanning tree

E條邊。選出V-1條。

[初始化]
答案設計成E個二進位數字,
0代表對應的邊不是MST。
1代表對應的邊是MST。

[fitness function]
是MST的權重,數值越大越好。

用人眼觀察,很直覺的就能在小區域找出最佳的子樹,那便是一套「特定組合」。

範例:travelling salesman problem

N個城市,C(N,2)條路。選出N條。

[初始化]
答案被迫設計成0到N-1的數字,代表一條路徑。

[fitness function]
是路徑的權重,數值越小越好。

[crossover]
哈哈哈,很難搞。

[mutation]
可以有好幾種方式:
1. 隨便交換兩個數字。
2. 隨便找一段數字,顛倒前後順序。
3. 隨便拿出一個數字,隨便插入到其他地方。

旅行推銷員問題也具有「特定組合」的性質,只不過它的答案並不適合分裂再結合,最好不要堅持使用基因演算法,自尋煩惱。

範例:video game

ant colony optimization(螞蟻演算法)

靈感來自於螞蟻覓食的過程,螞蟻發現食物後會沿途留下強烈的費洛蒙,驅使其他螞蟻沿著路線前進,不斷留下更多費洛蒙,吸引更多螞蟻;也有螞蟻會另闢新路,而找到更簡潔的路線。螞蟻演算法把答案比擬成螞蟻覓食的路徑,把好的答案不斷做局部調整,成為更好的答案。

螞蟻演算法有各式各樣的版本,這裡介紹其中一個經典版本ant colony system,主要是計算一條最短的覓食路徑。

附帶一提,這與螞蟻的真實行為完全無關。

1.
[初始化]
一開始先建立一個地圖,地圖上有P個地點。
有一些地點是食物,有一個地點是蟻窩。

地點與地點間皆有道路,
所有道路的費洛蒙都預設為1.0。

2.
[state transition rule]
一隻螞蟻從蟻窩出發。
每次踏上一個地點,螞蟻有兩種選擇,
  ◎探索:隨便走一條路。機率為q。
  ◎追蹤:走費洛蒙最強的道路。機率為1-q。

q是螞蟻選擇探索的機率,自行設定在0到1之間。

為了不讓螞蟻打轉繞圈,所以會讓螞蟻避開去過的地點。
在探索和追蹤前,要先過濾掉去過的地點。

所有食物都找到後就直接回蟻窩,
沒找到所有食物就不准回蟻窩。
總之就是要刻意弄出一條「嘗遍天下美食」的路線。

3.
[local updating rule]
螞蟻回巢之後,
剛剛走過的每一條道路,費洛蒙都會加強,
道路的費洛蒙會變成下列兩者相加後的數值,
  ◎自然蒸發,餘下:原本費洛蒙值,乘上1-ρ。
  ◎螞蟻路過,添加:道路起點所連接的道路當中,最大的那個費洛蒙值,乘上p。

ρ是費洛蒙蒸發的程度,自行設定在0到1之間。
通常會把參數設的很好,讓相加之後,費洛蒙會比原本增加一些,而不是變更少。

4.
有N隻螞蟻,依序做2. 3.這件事。

5.
[global updating rule]
N隻螞蟻回巢之後,
地圖上所有道路的費洛蒙都會蒸發一定比例。
  ◎自然蒸發,餘下:原本費洛蒙值,乘上1-α。

而目前的最佳路線,在蒸發之後,竟會神奇地額外增加一些費洛蒙。
  ◎最佳路線,添加:目前最佳解的數值,倒數後,再乘上α。

α是費洛蒙蒸發的程度,自行設定在0到1之間。

6.
2. 3. 4. 5.,重複執行R次。
途中不斷記錄最好的路線。
1. 初始化地圖與費洛蒙。
2. 以下動作執行R次:
  1. N隻螞蟻依序找路線,不是同時找路線:
    1. state transition rule,一隻螞蟻找一條路線。 
       此路線是由蟻窩出發,經過所有食物各一次,然後回到蟻窩。
    2. 記錄目前最佳路線。
    3. local updating rule,更新路線費洛蒙。
  2. global updating rule,更新所有道路費洛蒙。

演算法的過程,大致就是如此。細部的實作方式、參數的調校方式,隨人話虎爛。

螞蟻數量足夠豐富,可以避免進入局部極值。隨機探索用於增加答案豐富性,跳脫局部極值。

範例:travelling salesman problem

N個城市,C(N,2)條路。選出N條。

[初始化]
答案設計成0到N-1的數字,代表一條路徑。
地圖上每個地點都有食物。
地圖上可以任意挑一地點作為蟻窩。

當答案具有「特定排列」的性質,就適合採用螞蟻演算法。