第10章

重複の除去

(やさしい版) Pearls of Functional Algorithm Design(関数プログラミングによるアルゴリズム設計の真珠)

どんな問題?

リストの中から同じ要素をひとつだけ残して重複を消したい、というよくある操作を考える。ふつうの実装(Haskell 標準の nub)は「最初に出てきた順」で残すが、この章ではそれを少し変えて、結果の並びが辞書式順序でいちばん小さくなるように選び直したい。

たとえば 入力が “calculus” のとき、 どちらも「c, a, l, u, s の 5 文字を 1 回ずつ含む」点は同じ。並べ方だけが違う。
目標の速さ ふつうの nub は最悪 Θ(n2) 回の比較が要る。しかし今回の「辞書式最小」の版は、比較の下界がぐっと下がって Ω(n log n) しかいらない(Bloom & Wright 2003)。目標は Θ(n log n) のプログラムを組み立てること。ただし、道のりは長い。

はじめに

ふつうの nub のおさらい

Haskell 標準の nub は次のように書ける。

nub :: Eq a ⇒ [a] → [a] nub [ ] = [ ] nub (x : xs) = x : nub (xs \\ [x])
Dart // ふつうの nub:最初に出てきた順に重複を除く List<T> nubFirst<T>(List<T> xs) { if (xs.isEmpty) return <T>[]; final x = xs.first; final rest = xs.sublist(1).where((e) => e != x).toList(); return [x, ...nubFirst(rest)]; }

xs \\ ys は「xs から ys の各要素を取り除いた残り」を意味する。この定義では、要素は 最初に登場した順で並ぶ。

新しい nub にする

今回は次の型に変える。

nub :: Ord a ⇒ [a] → [a]

要求も「辞書式順序で最小の解を返せ」に変わる。Eq から Ord になっているのは、順序の比較が必須になったから。

位置か中身か 微妙な違いに注意。以前は「部分列としての位置」が辞書式最小だった。いまは「部分列そのもの」が辞書式最小になる。

最初の版

仕様(何をしたいかを式で書く)

nub = minimum · longest · filter nodups · subseqs

左から見て:

  1. subseqs:入力リストの部分列を全部作る。
  2. filter nodups:そのうち重複を含まないものだけ残す。
  3. longest:残ったものの中で最も長いものだけを取り出す。
  4. minimum:さらにその中で辞書式最小のものを選ぶ。

そこから導ける再帰版

この仕様から、次の再帰的定義が導ける(詳細な導出は練習問題として省略)。

nub [ ] = [ ] nub (x : xs) = if x ∉ xs then x : nub xs else (x : nub (xs \\ [x])) min (nub xs)
Dart // 辞書式最小の nub(素朴版・指数時間)。 // x が残りに出てこないなら選ぶしかない。出てくるなら「今使う」か「後で使う」の小さい方。 int lexCompare<T extends Comparable>(List<T> a, List<T> b) { final n = a.length < b.length ? a.length : b.length; for (var i = 0; i < n; i++) { final c = a[i].compareTo(b[i]); if (c != 0) return c; } return a.length.compareTo(b.length); } List<T> lexMin<T extends Comparable>(List<T> a, List<T> b) => lexCompare(a, b) <= 0 ? a : b; List<T> nubNaive<T extends Comparable>(List<T> xs) { if (xs.isEmpty) return <T>[]; final x = xs.first; final rest = xs.sublist(1); if (!rest.contains(x)) { return [x, ...nubNaive(rest)]; } final withX = [x, ...nubNaive(rest.where((e) => e != x).toList())]; final without = nubNaive(rest); return lexMin(withX, without); }

意味は素直で、「x が残りに出てこないなら選ぶしかない」「出てくるなら、x を今使うか後で使うか、二つの候補の小さい方を取る」というもの。

問題点 この再帰は場合分けごとに呼び出しが 2 倍に増える可能性があり、最悪で指数時間かかる。ここから Θ(n log n) までもっていくのが本章の仕事。

一般化

分割統治は素直には効かない

目標が Θ(n log n) なので、まず分割統治法を思いつく。すなわち次を満たす join を探したい。

nub (xs ++ ys) = join (nub xs) (nub ys) (10.1)

ところが、そんな join は存在しない。実際 (10.1) は次を同時に要求してしまう。

join “bca” “c” = nub “bcac” = “bac”
join “bca” “c” = nub “bcabc” = “abc”

左辺は同じなのに右辺が食い違う。したがって join は作れない。
この例からついでに、nubfoldl でも foldr でも書けないことも分かる。

hub という一般化を導入する

xxsyxsxy のときに nub(x:y:xs) を展開してみる。min の結合性と (x:) が min に分配することを使うと、

nub (x : y : xs) = x : y : nub (xs \\ [x, y]) min x : nub (xs \\ [x])) min y : nub (xs \\ [y]) min nub xs

ここで x < y なら第 3 項は第 2 項に負けるので落とせる:

nub (x : y : xs) = x : y : nub (xs \\ [x, y]) min x : nub (xs \\ [x])) min nub xs

逆に x > y なら第 1・2 項が負けるので落とせる:

nub (x : y : xs) = y : nub (xs \\ [y]) min nub xs

この二つの形は、「先頭にすでに決まった昇順の並び ws があり、そこに続けて残りをうまく選ぶ」という共通の骨格を示している。そこで hub を次のように定義する(以降 minimummin と略記)。

hub ws xs = min [is ++ nub (xs \\ is) | is ← inits ws] (10.2)

ここで ws厳密に増加するリスト、inits は接頭部(初期セグメント)を全部返す標準関数。上の例は次の形になる。

nub (x : y : xs) = if x < y then hub [x, y] xs else hub [y] xs
押さえておく事実

hub の帰納的定義を組み立てる

基底(リストが空)は簡単:

hub ws [ ] = {定義} min [is ++ nub ([ ] \\ is) | is ← inits ws] = {[ ] \\ is = [ ] かつ nub [ ] = [ ] より} min [is | is ← inits ws] = {[ ] は inits ws の中で辞書式最小} [ ]

よって hub ws [ ] = [ ]。

帰納段は次から始まる:

hub ws (x : xs) = min [is ++ nub ((x : xs) \\ is) | is ← inits ws] (10.3)

右辺を単純化するのに xws かどうかを知りたいので、wsx より小さい部分と残りに分ける:

(us, vs) = span (< x) ws

ws は増加順なので usvs も増加順。xws なら x = head vs、そうでなければ vs は空か x < head vs

さらに inits の次の性質が鍵となる:

inits (us ++ vs) = inits us ++ map (us ++) (inits+ vs) (10.4)

inits+ は空でない接頭部だけを返す。)これを (10.3) に代入して二つに分けると、hub ws (x : xs) = A min B が得られる。

A = min [is ++ nub ((x : xs) \\ is) | is ← inits us] (10.5)
B = min [us ++ is ++ nub ((x : xs) \\ (us ++ is)) | is ← inits+ vs] (10.6)

A の場合分け

A(xxs のとき)

A = {定義 (10.5)} min [is ++ nub ((x : xs) \\ is) | is ← inits us] = {x ∉ xs, x ∉ us より nub の再帰式から} min [is ++ [x] ++ nub (xs \\ is) | is ← inits us] = {is ∈ inits us に対し us < is ++ [x]} us ++ [x] ++ nub (xs \\ us) = {nub xs = hub [ ] xs} us ++ [x] ++ hub [ ] (xs \\ us)

A(xxs のとき)

A = {x ∈ xs, x ∉ us より nub の再帰式から} min [is ++ ([x] ++ nub (xs \\ (is ++ [x]))) min nub (xs \\ is)) | is ← inits us] = {min を内包表記の外へ} min [is ++ [x] ++ nub (xs \\ (is ++ [x])) | is ← inits us] min min [is ++ nub (xs \\ is) | is ← inits us] = {is ∈ inits us に対し us < is ++ [x]} (us ++ [x] ++ nub (xs \\ (us ++ [x]))) min min [is ++ nub (xs \\ is) | is ← inits us] = {inits (us ++ [x]) = inits us ++ [us ++ [x]]} min [is ++ nub (xs \\ is) | is ← inits (us ++ [x])] = {定義 (10.2)} hub (us ++ [x]) xs

まとめると、A は次と等しい:

if x ∈ xs then hub (us ++ [x]) xs else us ++ [x] ++ hub [ ] (xs \\ us)

B の展開

vs が空なら B は形式的な min [ ](比較対象なし)。そうでなければ:

B = {定義 (10.6)} min [us ++ is ++ nub ((x : xs) \\ (us ++ is)) | is ← inits+ (v : vs′)] = {inits+ (v : vs′) = map (v :) (inits vs′)} min [us ++ [v] ++ is ++ nub ((x : xs) \\ (us ++ [v] ++ is)) | is ← inits vs′] = {min · map (ys ++) = (ys ++) · min} us ++ [v] ++ min [is ++ nub ((x : xs) \\ (us ++ [v] ++ is)) | is ← inits vs′]

とくに Bus ++ [v] で始まる。

x が ws に入っていない場合

xws なら、vs は空か、または先頭 vx < v を満たす。どちらでも Aus ++ [x] の接頭部で始まり、us ++ [x] < us ++ [v] なので A < B。よって

x ∉ ws ⇒ hub ws (x : xs) = A

x が ws に入っている場合

このとき x = vB は次のように簡単化される:

us ++ [x] ++ min [is ++ nub (xs \\ (us ++ [x] ++ is)) | is ← inits vs′]

さらに xxsxxs で場合分け。

xxs のとき

B = {上より} us ++ [x] ++ min [is ++ nub (xs \\ (us ++ [x] ++ is)) | is ← inits vs′] = {xs \\ (us ++ [x] ++ is) = xs \\ (us ++ is) = (xs \\ us) \\ is} us ++ [x] ++ min [is ++ nub ((xs \\ us) \\ is) | is ← inits vs′] = {定義 (10.2)} us ++ [x] ++ hub vs′ (xs \\ us)
hub ws (x : xs) = {x ∉ xs を仮定した A と B の式より} (us ++ [x] ++ nub (xs \\ us)) min (us ++ [x] ++ hub vs′ (xs \\ us)) = {hub vs′ (xs \\ us) ≤ nub (xs \\ us)} us ++ [x] ++ hub vs′ (xs \\ us)

xxs のとき

hub ws (x : xs) = {x ∈ xs を仮定した A と B の式より} hub (us ++ [x]) xs min us ++ [x] ++ min [is ++ nub (xs \\ (us ++ [x] ++ is)) | is ← inits vs′] = {ws = us ++ [x] ++ vs′ と (10.4) より} min [is ++ nub (xs \\ is) | is ← inits ws] = {hub の定義} hub ws xs

ひとまずの完成形(図 10.1)

nub = hub [ ] hub ws [ ] = [ ] hub ws (x : xs) = case (x ∈ xs, x ∈ ws) of (False, False) → us ++ [x] ++ hub [ ] (xs \\ us) (False, True) → us ++ [x] ++ hub (tail vs) (xs \\ us) (True, False) → hub (us ++ [x]) xs (True, True) → hub ws xs where (us, vs) = span (< x) ws
Dart // ws は厳密増加、xs は残りのリスト。span (< x) ws を us / vs に分ける。 List<T> nubHubList<T extends Comparable>(List<T> xs) => hubList(<T>[], xs); List<T> hubList<T extends Comparable>(List<T> ws, List<T> xs) { if (xs.isEmpty) return <T>[]; final x = xs.first; final rest = xs.sublist(1); // (us, vs) = span (< x) ws var i = 0; while (i < ws.length && ws[i].compareTo(x) < 0) i++; final us = ws.sublist(0, i); final vs = ws.sublist(i); final xInXs = rest.contains(x); final xInWs = vs.isNotEmpty && vs.first.compareTo(x) == 0; final restMinusUs = rest.where((e) => !us.contains(e)).toList(); if (!xInXs && !xInWs) { return [...us, x, ...hubList(<T>[], restMinusUs)]; } else if (!xInXs && xInWs) { return [...us, x, ...hubList(vs.sublist(1), restMinusUs)]; } else if (xInXs && !xInWs) { return hubList([...us, x], rest); } else { return hubList(ws, rest); } }
図 10.1 nub の第二の定義
現状の計算量 所属判定・リスト差・span はいずれも線形時間。再帰呼び出しごとに線形かかるので、全体は 二乗時間。まだ目標には届かない。

集合の導入

Data.Set を使う

補助演算の計算量を線形から対数に落とすため、標準ライブラリ Data.Set を使う。使う演算は次のとおり。

empty :: Set a member :: Ord a ⇒ a → Set a → Bool insert :: Ord a ⇒ a → Set a → Set a split :: Ord a ⇒ a → Set a → (Set a, Set a) elems :: Ord a ⇒ Set a → [a]

前処理:各要素の「あと」に何があるかを持たせる

各要素 x について「それより後ろに出てくる要素の集合」を先に用意しておくと、所属判定が O(log n) で済む。

(x1, {x2, x3, … xn}), (x2, {x3, …, xn}), … (xn, { })
preprocess :: Ord a ⇒ [a] → [(a, Set a)] preprocess xs = zip xs (tail (scanr insert empty xs))
Dart import 'dart:collection'; // SplayTreeSet を Data.Set 相当として使う(対数時間で挿入・所属判定・分割ができる)。 // 各 xi に対して「後ろに出てくる要素の集合 {x_{i+1}, ..., x_n}」を対にして返す。 List<(T, SplayTreeSet<T>)> preprocess<T extends Comparable>(List<T> xs) { final n = xs.length; // scanr insert empty xs を末尾から作る final scan = List<SplayTreeSet<T>>.filled(n + 1, SplayTreeSet<T>()); scan[n] = SplayTreeSet<T>(); for (var i = n - 1; i >= 0; i--) { scan[i] = SplayTreeSet<T>.of(scan[i + 1])..add(xs[i]); } // zip xs (tail scan) return [for (var i = 0; i < n; i++) (xs[i], scan[i + 1])]; }

scanr insert empty は末尾から集合を積み上げていくので、全体で O(n log n) 時間で済む。

素朴に置き換えるとまだ遅い(図 10.2)

nub = hub empty · preprocess preprocess xs = zip xs (tail (scanr insert empty xs)) hub ws [ ] = [ ] hub ws ((x, xs) : xss) = case (member x xs, member x ws) of (False, False) → eus ++ [x] ++ hub empty yss (False, True) → eus ++ [x] ++ hub vs yss (True, False) → hub (insert x us) xss (True, True) → hub ws xss where (us, vs) = split x ws eus = elems us yss = [(x, xs) | (x, xs) ← xss, not (member x us)]
Dart // 集合を導入した版。所属判定・分割は O(log n) だが、yss の再構築が線形かかるので全体は O(n^2)。 (SplayTreeSet<T>, SplayTreeSet<T>) splitAt<T extends Comparable>( T x, SplayTreeSet<T> s) { final lo = SplayTreeSet<T>(); final hi = SplayTreeSet<T>(); for (final e in s) { if (e.compareTo(x) < 0) { lo.add(e); } else if (e.compareTo(x) > 0) { hi.add(e); } } return (lo, hi); } List<T> nubWithSet<T extends Comparable>(List<T> xs) => hubSet(SplayTreeSet<T>(), preprocess(xs)); List<T> hubSet<T extends Comparable>( SplayTreeSet<T> ws, List<(T, SplayTreeSet<T>)> xss) { if (xss.isEmpty) return <T>[]; final (x, xs) = xss.first; final rest = xss.sublist(1); final (us, vs) = splitAt(x, ws); final eus = us.toList(); // elems us(昇順) final xInXs = xs.contains(x); final xInWs = ws.contains(x); final yss = [ for (final e in rest) if (!us.contains(e.$1)) e, ]; if (!xInXs && !xInWs) { return [...eus, x, ...hubSet(SplayTreeSet<T>(), yss)]; } else if (!xInXs && xInWs) { return [...eus, x, ...hubSet(vs, yss)]; } else if (xInXs && !xInWs) { return hubSet(SplayTreeSet<T>.of(us)..add(x), rest); } else { return hubSet(ws, rest); } }
図 10.2 集合の導入
なぜまだ遅いのか yss の計算が呼び出しごとに Ω(n log m) かかり、合計で Ω(n2) になってしまう。例として入力
[1 .. n] ++ [j | j ← [1 .. n], j mod 3 ≠ 0]
を与えると、出力は [1 .. n]。3 の倍数のたびに ws から二要素ずつ吐き出されて([1, 2], [4, 5], [7, 8], …)、yss の総コストが n の二乗になる。

追加引数 ps で無駄を消す(図 10.3)

問題は「もう出力に使った要素」をリストから毎回フィルタし直している点にある。そこで hub に「これから飛ばすべき集合 ps」を追加引数として持たせ、hub′ を次のように定義する。

hub′ :: Set a → Set a → [(a, Set a)] → [a] hub′ ps ws xss = hub ws [(x, xs) | (x, xs) ← xss, x ∉ ps]
Dart // 「これから飛ばすべき集合 ps」を追加引数に持つ hub'。仕様レベルの定義。 List<T> hubPrimeSpec<T extends Comparable>( SplayTreeSet<T> ps, SplayTreeSet<T> ws, List<(T, SplayTreeSet<T>)> xss) => hubSet(ws, [ for (final e in xss) if (!ps.contains(e.$1)) e, ]);

これで、us の大きさを m とすると、除外リストの更新コスト qsO(m log n) で済む(以前の O(n log m) から改善)。この更新を行う us の大きさの合計は高々 n なので、全体は O(n log n) に収まる。

nub = hub′ empty empty · preprocess preprocess xs = zip xs (tail (scanr insert empty xs)) hub′ ps ws [ ] = [ ] hub′ ps ws ((x, xs) : xss) = if member x ps then hub′ ps ws xss else case (member x xs, member x ws) of (False, False) → eus ++ [x] ++ hub′ qs empty xss (False, True) → eus ++ [x] ++ hub′ qs vs xss (True, False) → hub′ ps (insert x us) xss (True, True) → hub′ ps ws xss where (us, vs) = split x ws eus = elems us qs = foldr insert ps eus
Dart // 最終版:O(n log n)。ps に「もう出力に使った要素」を溜め、xss からのフィルタを都度やらない。 List<T> nubFinal<T extends Comparable>(List<T> xs) => hubPrime( SplayTreeSet<T>(), SplayTreeSet<T>(), preprocess(xs)); List<T> hubPrime<T extends Comparable>(SplayTreeSet<T> ps, SplayTreeSet<T> ws, List<(T, SplayTreeSet<T>)> xss) { // 末尾再帰をループに展開して、深い再帰を避ける。 final out = <T>[]; while (xss.isNotEmpty) { final (x, xs) = xss.first; final rest = xss.sublist(1); if (ps.contains(x)) { xss = rest; continue; } final (us, vs) = splitAt(x, ws); final eus = us.toList(); final xInXs = xs.contains(x); final xInWs = ws.contains(x); if (!xInXs && !xInWs) { out.addAll(eus); out.add(x); // qs = foldr insert ps eus final qs = SplayTreeSet<T>.of(ps)..addAll(eus); ps = qs; ws = SplayTreeSet<T>(); xss = rest; } else if (!xInXs && xInWs) { out.addAll(eus); out.add(x); final qs = SplayTreeSet<T>.of(ps)..addAll(eus); ps = qs; ws = vs; xss = rest; } else if (xInXs && !xInWs) { ws = SplayTreeSet<T>.of(us)..add(x); xss = rest; } else { xss = rest; } } return out; }
図 10.3 最終版

結び

ここまで来るのに、かなりの手間を要した。しかも最終アルゴリズムは、美しいとも直感的とも言い難い。「これほど短く述べられる問題に、もっと単純な解があるのでは?」という疑念は残るが、少なくとも著者は見つけられていない。主計算は入り組み、場合分けにも悩まされた。

この章のまとめ とはいえ、戦略そのものはよくあるものである:
  1. まず問題を再帰的に定式化する。
  2. その再帰を、効率的に実装できる一般化された版に拡張する。
この二段構えは、多くの効率的アルゴリズムの導出でくり返し現れるパターンである。
Dart // 動作確認:"calculus" を辞書式最小の nub にすると "aclus"。 void main() { final input = 'calculus'.split(''); print(nubFirst(input).join()); // calus (ふつうの nub:最初に出てきた順) print(nubNaive(input).join()); // aclus (辞書式最小・素朴版) print(nubHubList(input).join()); // aclus (hub による第二の定義) print(nubWithSet(input).join()); // aclus (集合導入版、O(n^2)) print(nubFinal(input).join()); // aclus (最終版、O(n log n)) }

参考文献

Bloom, S. L. and Wright, R. S. (2003). Some lower bounds on comparison-based algorithms. Unpublished research paper. Department of Computer Science, Steven’s Institute of Technology, Hoboken, NJ, USA.