読者です 読者をやめる 読者になる 読者になる

MyEnigma

とあるエンジニアのブログです。#Robotics #Programing #C++ #Python #MATLAB #Vim #Mathematics #Book #Movie #Traveling #Mac #iPhone

ナップサック問題をPythonの最適化モデリングツールcvxpyで解いてみた

アルゴリズムまるごと学習ブック(日経BP Next ICT選書)

アルゴリズムまるごと学習ブック(日経BP Next ICT選書)

目次

 

はじめに

以前、

Pythonの最適化モデリングツールであるcvxpyを紹介しましたが、

myenigma.hatenablog.com

今回は、このcvxpyを使って、

有名な組み合わせ最適化問題である

ナップサック問題をcvxpyで解いてみました。

 

ナップサック問題とは?

ナップサック問題は、

下記の記事のように、有名な組み合わせ最適化問題の一つです。

ナップサック問題 - Wikipedia

2.5 ナップサック問題

 

言葉でナップサック問題を説明すると、

ある容量©のナップサックに、

複数の品物を入れる際に、

その品物のサイズ(S_i)と値段(V_i)を元に、

最も値段の総和が高くなるように、

ナップサックに入れる品物を選ぶ問題です。

 

たまにスーパーなどでやっている、

袋詰め放題サービスで、

何をどれだけ入れれば、コスパが高いか?

という問題を考えればわかりやすいかと思います。

 

上記の問題を、最適化の数式で表すと、

下記のようになります。

f:id:meison_amsl:20170410070655p:plain

 

今回、その品物をナップサックに入れるかどうかが重要なので、

それぞれの品物の種類に応じたバイナリ値(0: 入れない, 1: 入れる)

を格納したベクトルを最適化における変数としました。

上記の制約を守った上で、

コスト関数を最大化するベクトルxを求めるのが目的になります。

 

上記の最適化式を見ると分かる通り、

ナップサック問題は線形計画問題の一つになります。

myenigma.hatenablog.com

 

ちなみに、上記の問題のように、

変数が連続値ではなく、整数値やバイナリ値を取るような場合、

一般的なソルバーではなく、それらの問題を解くことができる

ECOS_BBのような特別なソルバーが必要になることが多いようです。

 

Pythonの最適化モデリングツールcvxpyでの解法

ナップサック問題は一般的に動的計画法を使って解くことが多いのですが、

今回は凸最適化の問題として、ソルバーを使って解いてみたいと思います。

 

今回は下記記事のナップサック問題の例題をcvxpyで解いてみました。

基本的には、前述の式をcvxpyで記述しただけです。

#
# A sample code to solve a knapsack_problem with cvxpy
#
# Author: Atsushi Sakai (@Atsushi_twi)
#

import cvxpy
import numpy as np

size = np.array([21, 11, 15, 9, 34, 25, 41, 52])
weight = np.array([22, 12, 16, 10, 35, 26, 42, 53])
capacity = 100

x = cvxpy.Bool(size.shape[0])
objective = cvxpy.Maximize(weight * x)
constraints = [capacity >= size * x]

prob = cvxpy.Problem(objective, constraints)
prob.solve(solver=cvxpy.ECOS_BB)
result = [round(ix[0, 0]) for ix in x.value]

print("status:", prob.status)
print("optimal value", prob.value)
print("result x:", result)

 

上記のコードを実行すると、

下記のようにちゃんとナップサック問題を解けていることがわかります。

(0,1,3,4,5の品物を選べば良い)

f:id:meison_amsl:20170410071901p:plain

 

ベンチマーク

上記のコードを少し改造して、

品物の種類を10000種類まで増やしてみました。

10000種類でも、自分のMacなら5秒ほどで解くことができました。

結構早いですね。

#
# A sample code to solve a knapsack_problem with cvxpy
#
# Author: Atsushi Sakai (@Atsushi_twi)
#

import cvxpy
import numpy as np
import random
import time

N_item = 10000

size = np.array([random.randint(1, 50) for i in range(N_item)])
weight = np.array([random.randint(1, 50) for i in range(N_item)])
capacity = 100

x = cvxpy.Bool(size.shape[0])
objective = cvxpy.Maximize(weight * x)
constraints = [capacity >= size * x]

start = time.time()
prob = cvxpy.Problem(objective, constraints)
prob.solve(solver=cvxpy.ECOS_BB)
elapsed_time = time.time() - start
print ("elapsed_time:{0}".format(elapsed_time) + "[sec]")
result = [round(ix[0, 0]) for ix in x.value]

print("status:", prob.status)
print("optimal value", prob.value)
print("result x:", result)

 

複数個、品物を選択出来る場合のナップサック問題の解法

前述の問題は、

与えられた品物を一個しか選べませんでしたが、

一種類の品物を複数個選べる問題も作れます。

 

そんな時は、Boolベクトルを正の整数のベクトルにすればOKです。

下記のようなコードで解くことできます。

#
# A sample code to solve a knapsack_problem with cvxpy
#
# Author: Atsushi Sakai (@Atsushi_twi)
#

import cvxpy
import numpy as np

size = np.array([21, 11, 15, 9, 34, 25, 41, 52])
weight = np.array([22, 12, 16, 10, 35, 26, 42, 53])
capacity = 100

x = cvxpy.Int(size.shape[0])
objective = cvxpy.Maximize(weight * x)
constraints = [capacity >= size * x]
constraints += [x >= 0]

prob = cvxpy.Problem(objective, constraints)
prob.solve(solver=cvxpy.ECOS_BB)
result = [round(ix[0, 0]) for ix in x.value]

print("status:", prob.status)
print("optimal value", prob.value)
print("size :", size * x.value)
print("result x:", result)

 

上記のコードを実行すると、

下記のように、4個目の商品を11個買うのがいいみたいです。

(正しいかは不明ですが..)

f:id:meison_amsl:20170410073818p:plain

 

GitHubリポジトリ

上記のコードはすべて下記のGitHubリポジトリでも公開しています。

github.com

 

参考資料

myenigma.hatenablog.com

myenigma.hatenablog.com

myenigma.hatenablog.com

myenigma.hatenablog.com

myenigma.hatenablog.com

アルゴリズムまるごと学習ブック(日経BP Next ICT選書)

アルゴリズムまるごと学習ブック(日経BP Next ICT選書)