ぱたへね

はてなダイアリーはrustの色分けができないのでこっちに来た

ナップサック問題の近似解を求める

Courseraの離散最適化コースを勉強していたところ、二週目の問題がこなせず一ヶ月以上停滞してしまいました。物の数(n)が30個くらいなら授業の方法で解けるのですが、40個で全然駄目になりました。練習問題は1000個とか10000個まであり厳密解でやるアプローチは見るからに無理、Webで調べても近似解の求め方はあまり見つからず苦労しました。

そんなとき、近似アルゴリズム―離散最適化問題への効果的アプローチ― という本がちょうど良いタイミングで発売され、そこに書いてあるアルゴリズムで上手く計算できたのでまとめます。

https://www.kyoritsu-pub.co.jp/bookdetail/9784320121775

Corseraの練習問題が元ネタです。そのまま回答になってしまうためコードは無しですいません。ちなみにRustで書きましたが、問題が大きくなると時間よりも先にメモリ不足で落ちてしまいました。そんなときも、近似解を求めることで使用メモリ量も制御できます。

問題

 入力:n組の2個の正整数 a_i, c_i(i=1,2,...,n) と正整数b

 タスク: \Sigma_{i\in X} a_i \leq b を満たすような部分集合 X \subseteq  {1,2,...n} のうち、和 \Sigma_{i \in X} ci が最大になる物を求める

ナップサックの用語にするとaiが大きさ、ciが価値、bがナップサックの大きさです。nの物品のなかから大きさの合計がb以下になる組み合わせで、価値の合計が最大になる組み合わせを探す問題です。

アルゴリズム(厳密解)

ナップサック問題における疑似多項式時間アルゴリズム(pseudo polynomial time algorithm)を使います。アルゴリズム自体は厳密解を求めます。最後に書くように、同じアルゴリズムで簡単に計算量を減らした近似解を求めることもできます。アルゴリズム時代はDP(動的計画法)の一種なのですが、いまいち正しい名前のような気がしません。

価値の最大値

 C = max \{ c_i : i = 1, 2, ..., n \}

とおくと   O(n^{2}C)  の計算時間で解くことができます。Cの値によって計算量が変わってきます。

今回はこの例でやってみます。

 a_i = \{4, 5, 8, 3, 2\}

 c_i = \{6, 10, 12, 4, 4\}

b = 11, n = 5です。

テーブルを作る

DP法と同じように二次元のテーブルを作ります。全ての物品の価値を足してPとします。

 P = \Sigma_{i=1}^{n} c_i = c(\{1, 2, ... n\})

配列をこのように作ります。

# 配列を初期化
a = [4, 5, 8, 3, 2];
c = [6, 10, 12, 4, 4];

テーブルの初期化

n×P+1の二次元配列を作る。 A[i][0]を0で初期化し、A[0][最初の物の価値]に最初の物の重さを入れる。それ以外は無限大とする。

配列の中を計算する

iとpで2重ループさせる。ループを疑似コードで書くとこう。

for(i=1;i<n;i++) {
  for(p=0;p<=P;p++) {
    処理
    }
}

処理の部分は、

  • p<c[i]ならば A[i][p] = A[i-1][p];
  • p >= ciならばA[i][p]は、①A[i-1][p]と②a[i]+A[i-1,p-c[i]]の小さい方とする。

良く分からないのでやってみましょう。 Pが、c = [6, 10, 12, 4, 4]の各要素の和なので36です。n=5なので、5x37の二次元配列Aを作ります。 iとpでループを回し、A[i][p]でアクセスします。

最初にi=0の行を初期化します。 A[0] = 0, A[c[0]] = a[0];で後の項はA[6]=4;になります。A[0]行のそれ以外の要素は無限大で初期化します。

p 0 1 2 3 4 5 6 7 8 ...
A[0][p] 0 4 ...

右側はp=36まで∞が続きます。

次にi=1の行を処理します。 c[i]=10なので、c[0]からc[9]までは上の列をコピーします。 c[10]は、①がA[0][10]で∞、②がa[1] + A[0][p-c[1]]=5+A[0][10-10] = 5+A[0][0] =5 になります。A[1][10]=min(①、②)なので、A[1][10]は5になります。同じように、A[1][16]が5+4で9になり、それ以外は∞です。

p 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 ...
A[0][p] 0 4 ...
A[1][p] 0 4 5 9 ...

これをi=4まで繰り返すとこうなります。

p 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36
A[0][p] 0 4
A[1][p] 0 4 5 9
A[2][p] 0 4 5 8 9 12 13 17
A[3][p] 0 3 4 5 8 8 9 12 12 13 16 17 20
A[4][p] 0 2 4 5 5 8 7 9 10 11 13 14 15 17 18 19 22

最大値探す

Aの値が決まると、次に一番下の行から、A[4][p] <= bとなる最大のpを探します。

p 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36
A[4][p] 0 2 4 5 5 8 7 9 10 11 13 14 15 17 18 19 22

bが11なので、A[4][20]、つまりp=20が最大値となります。この値が条件を満たす最大の価値qになります。

ナップサックの中身を決める

今度はq=pとして、iでn-1から1まで逆向きにループします。 疑似コードです。

for(i=n-1;i>=1;i--) {
  処理
}

処理の部分は、A[i][q]<A[i-1][q]ならば、そのiをナップサックの中に入れて、q=q-c[i]とします。ちょうど、表の一つ上と比較して、値が大きくなっていたらその物を追加し、qを価値分減らしていきます。最後、ループを出た後にq>0であれば、i=0の物もナップサックに加えます。厳密解を求めるときは、ループを出た後のqは、0かc[0]のどちらかです。

P=20より右側は無視してかまいません。

p 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
A[0][p] 0 4
A[1][p] 0 4 5 9
A[2][p] 0 4 5 8 9 12 12
A[3][p] 0 3 4 5 8 8 9 12 12
A[4][p] 0 2 4 5 5 8 7 9 10 11
  • i=4:A[4][20]の一つ上を見て11 < 12 なので、ナップサックに入れ、qを16にします。
  • i=3:A[3][16]の一つ上をみて、同じ値なのでナップサックに入れない。qはそのまま。
  • i=2:A[2][16]の一つ上をみて、同じ値なのでナップサックに入れない。qはそのまま。
  • i=1:A[1][16]の一つ上をみて、9 < ∞なのでナップサックにいれ、qを4にします。
  • 最後、q!=0なので、i=0もナップサックに入れます。

最終的に i={0,1,4}が厳密解になります。

アルゴリズム(近似解)

このアルゴリズムの計算量は [tex] O(n^{2}C) ] なので、Cが大きくなると遅くなります。例えば、C>nという状況では、n3以上の計算オーダーになります。

近似解の求め方は簡単で、物の重さ(c)とナップサックの容量bを、ある値で割って小さくする(Scaling)事で近似解を得ることができます。

a = [422, 543, 802, 322, 223];
b = 11234;

のような入力に対してaもbも100で割り最適解を求めると、それがそのまま近似解の組み合わせとなります。近似解の精度も保証されており、証明が気になる方は上の本を買って読んでみてください。