ここしばらくは「巡回セールスマン問題」について試行錯誤していました。きっかけは年末から開催されていたKaggleの「Traveling Santa 2018 - Prime Paths」コンペです（締め切り：2019年1月10日）。Kaggleでは年末恒例のイベントコンペのようで、一応賞金もでます（賞金総額：$25,000）。
ディープラーニングでは解けないようですが、面白そうなのでやってみることにしました（やってみれば得られることも多そうなので）。
巡回セールスマン問題（Traveling Salesman Problem/TSP）：
TSPとは、地図上にある複数の都市をセールスマンが一筆書きのように移動して、その合計移動距離が最小になる経路を求めるアルゴリズムです。地図上の都市が少なければ（約20箇所以下）、移動順番の組み合わせを総当たりで計算し最短経路を探し出すことができますが、都市数が増えるほど組み合わせの数は爆発的に増えてしまい、スーパーコンピュータを使っても数万年かかってしまうようです。すべての組み合わせを計算することができないため、近似計算、最適化などの方法によって少しでも正解に近づくように工夫するようです。TSPサイト。
今回のKaggleコンペルール：
基本はTSPにおける最短経路探索ですが、都市数は197769もあり、当然ながら総当たりで計算することが不可能。
都市名は0〜197768までのナンバリングされた数値であり、それぞれの都市には(X,Y)の座標値が割り当てられています。都市名0を開始点として都市間を移動し続け、最終的に巡回してまた0に戻るように都市の訪問順のリストをつくりあげて提出します。
さらに今回の特別ルールとして、都市から都市へ移動する際に、10ステップごとにペナルティが課されるときがあります。10ステップごとに一つ手前（9ステップ目）が素数の都市名であればペナルティなし、素数でない場合は9ステップ目から10ステップ目までの移動距離が1.1倍されてしまいます。要はできるだけ10ステップごとに9ステップ目が素数の都市名を通るように移動することが望ましいことになります。
実際どのくらいの都市数で、どのような経路になるかというと以下の画像のような感じです。
トナカイとソリの画像がそのまま都市の配置になっています。青い線が経路。これだけ都市数があるので、早々簡単には求められないということがわかります。
始めるにあたって：
まずTSPについては、一度は見たことがあるような問題ですが、詳しくは知らないのであくまでビギナーとして参加になります。調べていくとやはりTSPの研究者（数学者）がいるので、その人たちには到底敵いません。問題自体は至ってシンプルなのですが、解けそうで解けないという感じでしょうか。
KaggleのKernelsにはベースラインとなるコードが参加者によって掲示されているので、まずはそこをチェック。
TSP Solverライブラリ：
Kernelsを見るとTSP Solverというライブラリがいくつか存在し、それで計算すると総移動距離（スコア）が1533242前後になるようです。このスコアが少ないほど上位になるというわけです。多くの人は、このようなライブラリで数時間計算させベースとなる経路を得たのち、最適化アルゴリズムで少しづつ経路を調整してスコアを縮めていくようです。
アルゴリズム：
最適化アルゴリズムで少しずつスコアを縮めていくことはできますが、ある程度最適化してしまうと、それ以上スコアの変動はなくなってしまい手詰まりに陥ってしまいます。最初のTSPライブラリでかなりいいスコアを出しておかないときついという感じ。数日間連続計算させることでいい結果がでるときもありますが、今回はせいぜい一晩を上限としてプログラムを組んでみることにしてみました。組み合わせの問題でもあるので、組み合わせ数が多いほど有利になりますが、その分計算量も莫大になってしまいます。
とりあえずは、もう少しTSPの仕組みを理解しようとTSPライブラリを使わず自力でやってみることにしました。
欲張り法：
まずスタート地点から移動する際に、最も近い都市に移動し続けていけばいいのでは？と直感的に考えてしまいます。どうやら調べてみると、これは「欲張り法」というらしい。ということから、以下の関数が必要なので自前でつくることから始めてみました。
・二都市間の距離を測定する関数（単純距離）
・二都市間の移動経路のスコアを計算する関数（素数のペナルティも含め）
・ある都市に最も近い都市を求める関数（欲張り法用）
複素数で距離を求める：
(x1,y1)と(x2,y2)の二点間の距離を求めるには、
d = sqrt((x1-x2)**2 + (y1-y2)**2)
となりますが、pythonにおいては複素数がすぐに使えることから、y座標に1j（虚数）を掛けて各点を実部と虚部の一つの式にして計算できるようです。
p1 = x1 + y1*1j
p2 = x2 + y2*1j
d = abs(p1-p2)
都市が197768もあるので、順番に各都市との距離などを計算させると結構時間がかかります。できるだけ計算コストを抑えた工夫も必要そうです。上記の違いでどれだけ高速になっているかの検証はしていませんが、このような計算方法もあるというのが分かりました。
numbaで高速計算：
基本的にはnumpyで計算するためGPUを使えないことからどうしても計算が遅くなってしまいます。Kernelsを見ると、numbaを使えば高速計算や並行処理が可能になるらしく、慣れないながらも使ってみることにしました。numbaの最も簡単な使い方はdef関数前にデコレータとして@jitを追加させるだけです。
from numba import jit
@jit
def func(x):
y = some_calculation...
return y
このほかcudaを使う方法などもあるようですが、少々複雑なので今回は簡単な方法をできる範囲で使ってみることにしてみました。いずれにせよ、Kernelsからはこのような細かいテクニックが手に入るので結構勉強になります。
欲張り法の結果：
最短距離で繋いでいくので、そこそこいい結果がでるのではないかと期待していましたが、スコアは約1800000。TSPライブラリと比較すると約2割ほど多くなってしまいます。最初のうちは最も近い都市へ移動していくので効率的ですが、残りの都市が少なくなって行くに従って近くの都市も少なくなってしまうため、最短の都市であっても遠くの都市を選ばざるを得なくなってしまうからなのでしょう。
二方向の欲張り法：
欲張り法だと後半の都市選択が厳しくなってくることと、巡回ルートにしなければいけないことから、開始点から二方向に欲張り法でルートを繋げていくことにしてみました。一方が開始点からの順方向の移動で、もう一方が開始点までの逆方向の移動を最終的に中間地点でつなぐというやり方です。
しかしながらこの方法でもスコアは1800000程度で決定的といえるほど効果がありませんでした。単純に近い都市だけをつなぐだけだと進む方角の優先順位がないので、もう少し賢い方法が必要そうです。局所的に都市の分布密度が高い方へ優先的に移動するなどしたほうがよさそうですが、計算コストがかかりすぎるのでそこまではやらないことに。
ただ、あとからの最適化アルゴリズムでもスコアを縮めていくことは可能なので、とりあえず欲張り法で繋げておいてもいいのかもしれませんが、実際どうなのか（専門家ではないので分からない）？
このようにいろいろ試行錯誤しているとTSPの様々な問題が見えてきて、どうすればいいのか徐々に分かってきて興味も湧いてきます。
最適化アルゴリズム：
TSPライブラリであれ欲張り法であれ一旦経路をつくってから、それを調整していくアルゴリズムについてです。2-optや3-optあるいはk-optと言われるアルゴリズムがあるようです。
連続する二点入れ替え：
まず思いつくのが、任意の二点を入れ替えてみてスコアを計算し、もしスコア向上が得られれば採用するというアルゴリズム。このように何か最適化の処理をした後にスコアを評価して、結果がよければそれを採用するというやり方にいくつかあるようです。
最も簡単なのは連続する二点を入れ替える方法。
0 - 1 - 2 - 3 -(4)-(5)- 6 - 7 - 8 - 9 （初期状態）
4と5を入れ替える。
0 - 1 - 2 - 3 -(5)-(4)- 6 - 7 - 8 - 9 （変更後）
3〜6（4、5の前後）までの範囲のスコアを比較計算して向上すれば採用。
二点間反転：
また、遠く離れた二点間を入れ替える場合は、二点だけでなくその間の都市も含めて反転させてつなぎ直すという方法もあります。
0 - 1 - 2 -(3 - 4 - 5 - 6)- 7 - 8 - 9 （初期状態）
3から6を反転する。
0 - 1 - 2 -(6 - 5 - 4 - 3)- 7 - 8 - 9 （変更後）
2〜7（3〜6の前後）までのスコアを比較計算して向上すれば採用。
連続する三点入れ替え：
0 - 1 - 2 - 3 -(4)-(5)-(6)- 7 - 8 - 9 （初期状態）
(4, 5, 6)の三点の組み合わせすべて（6通り）を比較計算する。
0 - 1 - 2 - 3 -(4)-(5)-(6)- 7 - 8 - 9
0 - 1 - 2 - 3 -(4)-(6)-(5)- 7 - 8 - 9
0 - 1 - 2 - 3 -(5)-(4)-(6)- 7 - 8 - 9
0 - 1 - 2 - 3 -(5)-(6)-(4)- 7 - 8 - 9
0 - 1 - 2 - 3 -(6)-(4)-(5)- 7 - 8 - 9
0 - 1 - 2 - 3 -(6)-(5)-(4)- 7 - 8 - 9
3〜7までの範囲を比較計算し、最もスコアの少ない組み合わせを採用。
同様に四点、五点と増やしていけばスコアは向上するけれども、その分組み合わせのパターンも増大するので、ほどほどにしておかないと時間がかかりすぎる（せいぜい五点程度）。
組み合わせ用ライブラリ：
上記のような組み合わせのパターンを導き出すライブラリとしてsympy.utilities.iterables.multiset_permutationsがあります。これを使えば簡単に組み合わせパターンをリストとして取り出せます。
一点を任意の位置へ移動：
また別の方法としてシンプルなのが、一点を異なる場所へ入れ直す方法。
0 - 1 - 2 -(3)- 4 - 5 - 6 - 7 - 8 - 9 （初期状態）
3を1から順に挿入していき比較計算しスコア向上すれば更新する。
0 -(3)- 1 - 2 - 4 - 5 - 6 - 7 - 8 - 9
0 - 1 -(3)- 2 - 4 - 5 - 6 - 7 - 8 - 9
0 - 1 - 2 -(3)- 4 - 5 - 6 - 7 - 8 - 9
0 - 1 - 2 - 4 -(3)- 5 - 6 - 7 - 8 - 9
0 - 1 - 2 - 4 - 5 -(3)- 6 - 7 - 8 - 9
0 - 1 - 2 - 4 - 5 - 6 -(3)- 7 - 8 - 9
0 - 1 - 2 - 4 - 5 - 6 - 7 -(3)- 8 - 9
0 - 1 - 2 - 4 - 5 - 6 - 7 - 8 -(3)- 9
0 - 1 - 2 - 4 - 5 - 6 - 7 - 8 - 9 -(3)
アルゴリズムとしてはシンプルですが、これを開始点0を除いた全ての都市（197768箇所）において順に比較計算させていくと197768*197767回計算させることになり、かなり時間がかかってしまうので、移動したい点に最も近い点を上位五点まで選択しておき、それらのうちでスコア向上が見られるときにだけ採用とするなどの工夫をしました。
素数用ライブラリ：
この他、このコンペ特有のルールとなる素数に関しては、sympy.isprimeをつかうことで、その数が素数かどうかをすぐに判定できます。
最適化の結果：
上記の方法だけでもある程度最適化はできましたが、やはり最初のTSPライブラリによる結果がいいほうが最終的なスコアがよくなってしまいます。自前の欲張り法による初期経路の場合は、一度の最適化でもかなりスコアを向上させることができましたが（もともと改善の余地が多いため）、TSPライブラリによる初期経路の場合は、すでにかなり改善されているためか1%も向上しません。数時間回してやっと１箇所（0.0001%程度）改善するくらいでしょうか。
もっと強力なコンピュータか計算の高速化の工夫も必要そうなので、単なるTSPだけの問題というわけでもなくなってきて奥が深いという感じ。今回はPythonを使いましたが、C++でやったほうがよさそうです。
この他にもいくつか最適化するアルゴリズムを考えてみましたが、計算コストがかかりすぎることから諦めたものもあります（一通り計算させると１ヶ月以上かかってしまうなど）。
順位的には大したことはありませんでしたが、今回もまたTSPという未経験の分野にチャレンジしたことでいろいろと収穫がありました。やはりKaggleに参加するとモチベーションもあがり勉強にもなります。もうコンペ自体は終わってしまいましたが、引き続きTSPの最適化アルゴリズムを考えてみたくなります。
TSP：その２へ続く。
関連：Traveling Salesman Problem:巡回セールスマン問題について（まとめ）
ディープラーニングでは解けないようですが、面白そうなのでやってみることにしました（やってみれば得られることも多そうなので）。
巡回セールスマン問題（Traveling Salesman Problem/TSP）：
TSPとは、地図上にある複数の都市をセールスマンが一筆書きのように移動して、その合計移動距離が最小になる経路を求めるアルゴリズムです。地図上の都市が少なければ（約20箇所以下）、移動順番の組み合わせを総当たりで計算し最短経路を探し出すことができますが、都市数が増えるほど組み合わせの数は爆発的に増えてしまい、スーパーコンピュータを使っても数万年かかってしまうようです。すべての組み合わせを計算することができないため、近似計算、最適化などの方法によって少しでも正解に近づくように工夫するようです。TSPサイト。
今回のKaggleコンペルール：
基本はTSPにおける最短経路探索ですが、都市数は197769もあり、当然ながら総当たりで計算することが不可能。
都市名は0〜197768までのナンバリングされた数値であり、それぞれの都市には(X,Y)の座標値が割り当てられています。都市名0を開始点として都市間を移動し続け、最終的に巡回してまた0に戻るように都市の訪問順のリストをつくりあげて提出します。
さらに今回の特別ルールとして、都市から都市へ移動する際に、10ステップごとにペナルティが課されるときがあります。10ステップごとに一つ手前（9ステップ目）が素数の都市名であればペナルティなし、素数でない場合は9ステップ目から10ステップ目までの移動距離が1.1倍されてしまいます。要はできるだけ10ステップごとに9ステップ目が素数の都市名を通るように移動することが望ましいことになります。
実際どのくらいの都市数で、どのような経路になるかというと以下の画像のような感じです。
始めるにあたって：
まずTSPについては、一度は見たことがあるような問題ですが、詳しくは知らないのであくまでビギナーとして参加になります。調べていくとやはりTSPの研究者（数学者）がいるので、その人たちには到底敵いません。問題自体は至ってシンプルなのですが、解けそうで解けないという感じでしょうか。
KaggleのKernelsにはベースラインとなるコードが参加者によって掲示されているので、まずはそこをチェック。
TSP Solverライブラリ：
Kernelsを見るとTSP Solverというライブラリがいくつか存在し、それで計算すると総移動距離（スコア）が1533242前後になるようです。このスコアが少ないほど上位になるというわけです。多くの人は、このようなライブラリで数時間計算させベースとなる経路を得たのち、最適化アルゴリズムで少しづつ経路を調整してスコアを縮めていくようです。
アルゴリズム：
最適化アルゴリズムで少しずつスコアを縮めていくことはできますが、ある程度最適化してしまうと、それ以上スコアの変動はなくなってしまい手詰まりに陥ってしまいます。最初のTSPライブラリでかなりいいスコアを出しておかないときついという感じ。数日間連続計算させることでいい結果がでるときもありますが、今回はせいぜい一晩を上限としてプログラムを組んでみることにしてみました。組み合わせの問題でもあるので、組み合わせ数が多いほど有利になりますが、その分計算量も莫大になってしまいます。
とりあえずは、もう少しTSPの仕組みを理解しようとTSPライブラリを使わず自力でやってみることにしました。
欲張り法：
まずスタート地点から移動する際に、最も近い都市に移動し続けていけばいいのでは？と直感的に考えてしまいます。どうやら調べてみると、これは「欲張り法」というらしい。ということから、以下の関数が必要なので自前でつくることから始めてみました。
・二都市間の距離を測定する関数（単純距離）
・二都市間の移動経路のスコアを計算する関数（素数のペナルティも含め）
・ある都市に最も近い都市を求める関数（欲張り法用）
複素数で距離を求める：
(x1,y1)と(x2,y2)の二点間の距離を求めるには、
d = sqrt((x1-x2)**2 + (y1-y2)**2)
となりますが、pythonにおいては複素数がすぐに使えることから、y座標に1j（虚数）を掛けて各点を実部と虚部の一つの式にして計算できるようです。
p1 = x1 + y1*1j
p2 = x2 + y2*1j
d = abs(p1-p2)
都市が197768もあるので、順番に各都市との距離などを計算させると結構時間がかかります。できるだけ計算コストを抑えた工夫も必要そうです。上記の違いでどれだけ高速になっているかの検証はしていませんが、このような計算方法もあるというのが分かりました。
numbaで高速計算：
基本的にはnumpyで計算するためGPUを使えないことからどうしても計算が遅くなってしまいます。Kernelsを見ると、numbaを使えば高速計算や並行処理が可能になるらしく、慣れないながらも使ってみることにしました。numbaの最も簡単な使い方はdef関数前にデコレータとして@jitを追加させるだけです。
from numba import jit
@jit
def func(x):
y = some_calculation...
return y
このほかcudaを使う方法などもあるようですが、少々複雑なので今回は簡単な方法をできる範囲で使ってみることにしてみました。いずれにせよ、Kernelsからはこのような細かいテクニックが手に入るので結構勉強になります。
欲張り法の結果：
最短距離で繋いでいくので、そこそこいい結果がでるのではないかと期待していましたが、スコアは約1800000。TSPライブラリと比較すると約2割ほど多くなってしまいます。最初のうちは最も近い都市へ移動していくので効率的ですが、残りの都市が少なくなって行くに従って近くの都市も少なくなってしまうため、最短の都市であっても遠くの都市を選ばざるを得なくなってしまうからなのでしょう。
二方向の欲張り法：
欲張り法だと後半の都市選択が厳しくなってくることと、巡回ルートにしなければいけないことから、開始点から二方向に欲張り法でルートを繋げていくことにしてみました。一方が開始点からの順方向の移動で、もう一方が開始点までの逆方向の移動を最終的に中間地点でつなぐというやり方です。
しかしながらこの方法でもスコアは1800000程度で決定的といえるほど効果がありませんでした。単純に近い都市だけをつなぐだけだと進む方角の優先順位がないので、もう少し賢い方法が必要そうです。局所的に都市の分布密度が高い方へ優先的に移動するなどしたほうがよさそうですが、計算コストがかかりすぎるのでそこまではやらないことに。
ただ、あとからの最適化アルゴリズムでもスコアを縮めていくことは可能なので、とりあえず欲張り法で繋げておいてもいいのかもしれませんが、実際どうなのか（専門家ではないので分からない）？
このようにいろいろ試行錯誤しているとTSPの様々な問題が見えてきて、どうすればいいのか徐々に分かってきて興味も湧いてきます。
最適化アルゴリズム：
TSPライブラリであれ欲張り法であれ一旦経路をつくってから、それを調整していくアルゴリズムについてです。2-optや3-optあるいはk-optと言われるアルゴリズムがあるようです。
連続する二点入れ替え：
まず思いつくのが、任意の二点を入れ替えてみてスコアを計算し、もしスコア向上が得られれば採用するというアルゴリズム。このように何か最適化の処理をした後にスコアを評価して、結果がよければそれを採用するというやり方にいくつかあるようです。
最も簡単なのは連続する二点を入れ替える方法。
0 - 1 - 2 - 3 -(4)-(5)- 6 - 7 - 8 - 9 （初期状態）
4と5を入れ替える。
0 - 1 - 2 - 3 -(5)-(4)- 6 - 7 - 8 - 9 （変更後）
3〜6（4、5の前後）までの範囲のスコアを比較計算して向上すれば採用。
二点間反転：
また、遠く離れた二点間を入れ替える場合は、二点だけでなくその間の都市も含めて反転させてつなぎ直すという方法もあります。
0 - 1 - 2 -(3 - 4 - 5 - 6)- 7 - 8 - 9 （初期状態）
3から6を反転する。
0 - 1 - 2 -(6 - 5 - 4 - 3)- 7 - 8 - 9 （変更後）
2〜7（3〜6の前後）までのスコアを比較計算して向上すれば採用。
連続する三点入れ替え：
0 - 1 - 2 - 3 -(4)-(5)-(6)- 7 - 8 - 9 （初期状態）
(4, 5, 6)の三点の組み合わせすべて（6通り）を比較計算する。
0 - 1 - 2 - 3 -(4)-(5)-(6)- 7 - 8 - 9
0 - 1 - 2 - 3 -(4)-(6)-(5)- 7 - 8 - 9
0 - 1 - 2 - 3 -(5)-(4)-(6)- 7 - 8 - 9
0 - 1 - 2 - 3 -(5)-(6)-(4)- 7 - 8 - 9
0 - 1 - 2 - 3 -(6)-(4)-(5)- 7 - 8 - 9
0 - 1 - 2 - 3 -(6)-(5)-(4)- 7 - 8 - 9
3〜7までの範囲を比較計算し、最もスコアの少ない組み合わせを採用。
同様に四点、五点と増やしていけばスコアは向上するけれども、その分組み合わせのパターンも増大するので、ほどほどにしておかないと時間がかかりすぎる（せいぜい五点程度）。
組み合わせ用ライブラリ：
上記のような組み合わせのパターンを導き出すライブラリとしてsympy.utilities.iterables.multiset_permutationsがあります。これを使えば簡単に組み合わせパターンをリストとして取り出せます。
一点を任意の位置へ移動：
また別の方法としてシンプルなのが、一点を異なる場所へ入れ直す方法。
0 - 1 - 2 -(3)- 4 - 5 - 6 - 7 - 8 - 9 （初期状態）
3を1から順に挿入していき比較計算しスコア向上すれば更新する。
0 -(3)- 1 - 2 - 4 - 5 - 6 - 7 - 8 - 9
0 - 1 -(3)- 2 - 4 - 5 - 6 - 7 - 8 - 9
0 - 1 - 2 -(3)- 4 - 5 - 6 - 7 - 8 - 9
0 - 1 - 2 - 4 -(3)- 5 - 6 - 7 - 8 - 9
0 - 1 - 2 - 4 - 5 -(3)- 6 - 7 - 8 - 9
0 - 1 - 2 - 4 - 5 - 6 -(3)- 7 - 8 - 9
0 - 1 - 2 - 4 - 5 - 6 - 7 -(3)- 8 - 9
0 - 1 - 2 - 4 - 5 - 6 - 7 - 8 -(3)- 9
0 - 1 - 2 - 4 - 5 - 6 - 7 - 8 - 9 -(3)
アルゴリズムとしてはシンプルですが、これを開始点0を除いた全ての都市（197768箇所）において順に比較計算させていくと197768*197767回計算させることになり、かなり時間がかかってしまうので、移動したい点に最も近い点を上位五点まで選択しておき、それらのうちでスコア向上が見られるときにだけ採用とするなどの工夫をしました。
素数用ライブラリ：
この他、このコンペ特有のルールとなる素数に関しては、sympy.isprimeをつかうことで、その数が素数かどうかをすぐに判定できます。
最適化の結果：
上記の方法だけでもある程度最適化はできましたが、やはり最初のTSPライブラリによる結果がいいほうが最終的なスコアがよくなってしまいます。自前の欲張り法による初期経路の場合は、一度の最適化でもかなりスコアを向上させることができましたが（もともと改善の余地が多いため）、TSPライブラリによる初期経路の場合は、すでにかなり改善されているためか1%も向上しません。数時間回してやっと１箇所（0.0001%程度）改善するくらいでしょうか。
もっと強力なコンピュータか計算の高速化の工夫も必要そうなので、単なるTSPだけの問題というわけでもなくなってきて奥が深いという感じ。今回はPythonを使いましたが、C++でやったほうがよさそうです。
この他にもいくつか最適化するアルゴリズムを考えてみましたが、計算コストがかかりすぎることから諦めたものもあります（一通り計算させると１ヶ月以上かかってしまうなど）。
順位的には大したことはありませんでしたが、今回もまたTSPという未経験の分野にチャレンジしたことでいろいろと収穫がありました。やはりKaggleに参加するとモチベーションもあがり勉強にもなります。もうコンペ自体は終わってしまいましたが、引き続きTSPの最適化アルゴリズムを考えてみたくなります。
TSP：その２へ続く。
関連：Traveling Salesman Problem:巡回セールスマン問題について（まとめ）
門脇 大輔 阪田 隆司 保坂 桂佑 平松 雄司
技術評論社
売り上げランキング: 363
技術評論社
売り上げランキング: 363
