焼きなまし法
Pythonで焼きなまし法 (Simulated Annealing) をやってみる.
焼きなまし法は最適化問題を解くための乱択アルゴリズムの一つ.
状態をランダムに近傍の状態に変化させ,コストが小さくなる場合に状態を更新する.
ただし,状態遷移によってコストが増加する場合でも以下の確率 $ pで遷移を行う:
$ \Delta = f(x') - f(x)
$ p = \exp(-\Delta / T)
$ x, x'はそれぞれ現在の状態,次の状態, $ f はコスト関数,$ Tは温度パラメータ.
状態遷移を繰り返しながら徐々に温度$ Tを下げていくことで解を収束させる.
パラメータについて
焼きなまし法では,状態近傍,遷移確率,焼きなましスケジュールをパラメータとして設定する必要がある.
最も重要なのは状態近傍であり,コストが近くなるような近傍を設定することが望ましい.
遷移確率については上記に示した定義が使われる場合が多いが,数学的根拠があるわけではない.
焼きなましスケジュールについては初期状態において上りと下りの確率がほぼ等しくなる程度に設定する必要がある.
実装
ここでは以下のような単純なソルバーを用いる.
transit と cost メソッドを実装すると焼きなまし方が利用できる.
transitは次の状態を返すメソッドであり,costは状態からコストを計算するメソッドである.
code: annealer.py
import numpy as np
class SimulatedAnnealer:
def __init__(self, temperature: float = 10000, cool: float = 0.99) -> None:
self._temperature = temperature
self._cool = cool
def transit(self, state: np.ndarray) -> np.ndarray:
raise NotImplementedError
def cost(self, state: np.ndarray) -> float:
raise NotImplementedError
def solve(self, state: np.ndarray) -> np.ndarray:
curr_state = state
temperature = self._temperature
while temperature > 0.0001:
next_state = self.transit(curr_state)
curr_cost = self.cost(curr_state)
next_cost = self.cost(next_state)
transit_prob = np.exp(-abs(next_cost - curr_cost) / temperature)
if next_cost < curr_cost or np.random.uniform() < transit_prob:
curr_state = next_state
temperature *= self._cool
return curr_state
多項式の最小化
連続な実数値を扱うコスト関数$ f: \mathbb{R} \rightarrow \mathbb{R} を最小化する問題を解く.
まずはtransit と cost メソッドを実装する.
transitメソッドでは正または負の方向に少しだけ値を変化させる.
cost関数は引数として受け取れるようにした.
code: polynomial.py
class NumericalOptimizer(SimulatedAnnealer):
def __init__(
self,
stride: float = 1.0,
**kwargs,
) -> None:
super().__init__(**kwargs)
self._stride = stride
self._cost = cost
def transit(self, state: np.ndarray) -> np.ndarray:
fluct = np.random.uniform(size=state.shape)
delta = (fluct - 0.5) * self._strid
return state + delta
def cost(self, state: np.ndarray) -> float:
return self._cost(state)
ここでは以下のように $ x = 0.8 が大域的最適解となるような多項式のコスト関数を最小化する.
https://gyazo.com/9d231fafbc49098b073862f49ba1154f
code:cost.py
def cost(x):
a, b = 0.37, 0.8
return 3 * x**4 - 4 * (a + b) * x**3 + 6 * a * b * x**2
100回焼きなましを行い,得られた解をプロットする.
code:solve_polynomial.py
# optimize
num_trial = 100
results = np.zeros(num_trial)
for i in tqdm(range(num_trial)):
init_state = np.random.uniform(size=())
optimizer = NumericalOptimizer(cost=cost, stride=1)
resultsi = float(optimizer.solve(init_state)) # report summary
is_success = np.isclose(results, 0.8, atol=0.05)
print(f"success ratio: {is_success.mean()} "
f"({is_success.sum()} / {len(is_success)})"
# plot results
x = np.linspace(min(results) - 0.1, max(results) + 0.1, 1000)
plt.plot(x, cost(x))
color="orange",
label="success")
color="blue",
label="fail")
plt.legend()
plt.show()
https://gyazo.com/bbaf42e99e7a4c0db2be8aa09bb15864
焼きなましの結果,大域的最適解の0.8と局所的最適解の0.0の周辺に解が集まっていることがわかる.
巡回セールスマン問題
焼きなまし法を使って巡回セールスマン問題を解く.
巡回セールスマン問題では複数の都市を1回ずつ通って戻ってくる際に最短となる都市の順番を求める.
まずはtransit と cost メソッドを実装する.
transitではstateとして都市の順序を受け取り,このうち2つを交換して返す.
costではstateと距離行列に基づいて都市を巡回した場合の総距離を計算する.
code: tsp.py
class TSPOptimizer(SimulatedAnnealer):
def __init__(self, distance_matrix: np.ndarray, **kwargs) -> None:
super().__init__(**kwargs)
self._distance_matrix = distance_matrix
def transit(self, state: np.ndarray) -> np.ndarray:
state = state.copy()
a, b = np.random.choice(range(1, len(state) - 1), 2, replace=False)
statea, stateb = stateb, statea return state
def cost(self, state: np.ndarray) -> float:
total_distance = 0.0
num_nodes = len(state)
for indice in range(num_nodes):
total_distance += distance
return total_distance
最適化&可視化のための処理は以下:
code: solve_tsp.py
num_nodes = 10
nodes = build_nodes(num_nodes) # return matrix with shape (num_nodes, 2)
# compute distance matrix
distance_matrix = np.zeros((num_nodes, num_nodes))
for i in range(num_nodes):
for j in range(i, num_nodes):
distance = np.linalg.norm(nodesi - nodesj) distance_matrixij = distance distance_matrixji = distance # optimize
init_state = np.random.permutation(range(num_nodes))
optimizer = TSPOptimizer(distance_matrix)
result = optimizer.solve(init_state)
# plot result
plt.plot(
"--",
color="blue",
alpha=0.2,
label="initial state",
)
plt.plot(
color="orange",
label="result",
)
plt.scatter(
color="orange",
label="nodes",
)
plt.legend()
plt.show()
円周上に等間隔に都市を配置して焼きなまし法による最適化を行った:
code: build_circle_ndoes.py
def build_nodes(num_nodes: int) -> np.ndarray:
nodes = np.zeros((num_nodes, 2))
theta = np.linspace(0.0, 2.0 * np.pi, num_nodes)
nodes:, 0, nodes:, 1 = np.cos(theta), np.sin(theta) return nodes
https://gyazo.com/0446a55116514fa1c692d3937690c748
点は都市,破線は初期状態,実線は最適化結果を表す.
焼きなましの結果,最適解が得られていることがわかる.
次に都市をランダムに配置して焼きなましを行った:
code: build_random_nodes.py
def build_nodes(num_nodes: int) -> np.ndarray:
nodes = np.random.uniform(size=(num_nodes, 2))
return nodes
https://gyazo.com/8663e47c1c491f5e8d7665013df6bd00
ランダムに配置した場合でも,最適に近い結果が得られている.
参考