組合せ最適化で4色問題を解く

# <i class='fa fa-image' /> 4色問題を解く
[組合せ最適化](http://qiita.com/Tsutomu-KKE@github/items/bfbf4c185ed7004b5721)を使うと4色問題も解けます。
ここでは、画像から自動的に問題を読み取ってみましょう。問題は、[ニコリ様](http://www.nikoli.com/ja/take_a_break/four_color_problem/)からいただきました。


# <i class='fa fa-image' /> 定式化
複数のエリアがあり、隣接したエリア同士は、異なる色で塗ることとします。

<table>
<tr><td>$\mbox{variables}$</td><td>$v_{ij} \in \{0, 1\} ~ \forall i, j$</td><td>エリアiが色jかどうか (1)</td></tr>
<tr><td>$\mbox{subject to}$</td><td>$\sum_j{v_{ij}} = 1 ~ \forall i$</td><td>色を1つ選ぶ (2)</td></tr>
<tr><td></td><td colspan="2">隣り合うエリアは異なる色にする (3)</td></tr>
</table>


# <i class='fa fa-image' /> Pythonで解く
画像ファイルを読込みます。

```py3:python
import networkx as nx
from urllib import request
from PIL import Image, ImageDraw
from collections import Counter
from itertools import product
from random import seed, shuffle
from pulp import *

# 画像ファイルの読込み
with request.urlopen('https://dl.dropboxusercontent.com/u/35689878/tmp/4c.png') as fd:
    im = Image.open(fd)
```

![4c.png](https://qiita-image-store.s3.amazonaws.com/0/13955/8df50989-747d-cc7f-a28d-e2d6691a6fa0.png)


- ターゲットとなる代表色を頻度の多い色として抽出します。代表色で塗りつぶされたものがエリアです。
- 結果の色でエリアを塗るために、一旦、エリアをわかりやすい色(赤が0、緑が1、青がエリア番号)で塗ることにします。
- 最初に RGB=(0,1,?)の色のG(緑)を0にします。
- 代表色のエリアをRGB=(0,1,通し番号)で塗りつぶします。
- 境界を少し広げます
  - ランダムにポイントを選び、その色がエリアで1ピクセル隣りがどのエリアでもなければそのエリアにします。

```py3:python
# 代表色(最も使用頻度の多い色)を抽出
cc = sorted([(v, k) for k, v in Counter(im.getdata()).items()])[-1][1]

# RGB=(0,1,?)の色をなくす
for y, x in product(range(im.height), range(im.width)):
    R, G, B = im.getpixel((x, y))[:3]
    if (R, G) == (0, 1):
        im.putpixel(0, 0, B)

# 代表色のエリアをRGB=(0,1,通し番号)で塗りつぶす
n = 0
for y, x in product(range(im.height), range(im.width)):
    if im.getpixel((x, y)) != cc:
        continue
    ImageDraw.floodfill(im, (x, y), (0, 1, n))
    n += 1

# 境界を少し広げる
seed(1)
dd = [(-1, 0), (0, -1), (0, 1), (1, 0)]
for h in range(1):
    l = list(product(range(1, im.height-1), range(1, im.width-1)))
    shuffle(l)
    for y, x in l:
        c = im.getpixel((x, y))
        if c[:2] == (0, 1):
            for i, j in dd:
                if im.getpixel((x+i, y+j))[:2] != (0, 1):
                    im.putpixel((x+i, y+j), c)
```

エリアを点とし、ピクセルが隣り合うエリアの間を辺とするグラフを作成します。

```py3:python
# 隣り合うかどうかをグラフで表す
g = nx.Graph()
for y, x in product(range(im.height-1), range(im.width-1)):
    c1 = im.getpixel((x, y))
    if c1[:2] != (0, 1):
        continue
    c2 = im.getpixel((x+1, y))
    c3 = im.getpixel((x, y+1))
    if c2[:2] == (0, 1) and c1[2] != c2[2]:
        g.add_edge(c1[2], c2[2])
    if c3[:2] == (0, 1) and c1[2] != c3[2]:
        g.add_edge(c1[2], c3[2])
```

定式化して解きます。

```py3:python
# 4色問題を解く
r4 = range(4)
m = LpProblem() # 数理モデル
# エリアiを色jにするかどうか (1)
v = [[LpVariable('v%d_%d'%(i, j), cat=LpBinary) for j in r4] for i in g.nodes()]
for i in g.nodes():
    m += lpSum(v[i]) == 1 # (2)
for i, j in g.edges():
    for k in r4:
        m += v[i][k] + v[j][k] <= 1 # (3)
m.solve()
co = [(97, 132, 219), (228, 128, 109), (255, 241, 164), (121, 201, 164)] # 4色
rr = [int(value(lpDot(r4, i))) for i in v] # 結果
for y, x in product(range(im.height-1), range(im.width-1)):
    c = im.getpixel((x, y))
    if c[:2] == (0, 1): # エリアならば、結果で塗る
        ImageDraw.floodfill(im, (x, y), co[rr[c[2]]])
im.save('result.png')
```

![result.png](https://qiita-image-store.s3.amazonaws.com/0/13955/8d95f45c-2d32-e080-d75f-63d70436fd86.png)


実際に試してみたらOKでした。

![cong.png](https://qiita-image-store.s3.amazonaws.com/0/13955/946741a6-88ea-fe03-902b-f2097a8fb38d.png)

以上