初夏のワークショップ

June 15, 2016

おねえさん問題の乱択近似

柴田友樹1, 山内由紀子2, 来嶋秀治2, 山下雅史2
九州大学 工学部電気情報工学科 (1:研究当時)
大学院システム情報科学研究院 (2)[参照] 柴田友樹，山内由紀子, 来嶋秀治, 山下雅史, 格子上の経路数え上げの乱択近似--24時間で解くおねえさん問題, 2015年度 冬のLAシンポ
ジウム, 京都大学数理解析研究所, 2016年1月28日.
 柴田 友樹, 格子上の経路数え上げの乱択近似, 九州大学工学部電気情報工学科卒業論文, 2016年3月.

2おねえさん問題
格子グラフ上で左上(𝑠𝑠)から右下(𝑡𝑡)を結ぶ，𝑠𝑠-𝑡𝑡単純経路の総数を求めよ．

画像引用 ： 日本科学未来館, 『フカシギの数え方』 おねえさんといっしょ！
みんなで数えてみよう！, YouTube, 2012/9/10公開,
https://www.youtube.com/watch?v=Q4gTV4r0zRs*H. Iwashita, Y. Nakazawa, J. Kuwahara, T. Uno, S. Minato,Efficient computation of the number of paths in a grid graphwith minimal perfect hash functions, Hokkaido UniversityTCS Technical Report, TCS-TR-A-13-64, 2013.

26 17369931586279272931175440421236498900372229588288140604663703720910342413276134762789218193498006107082296223143380491348290026721931129627708738890853908108906396

25 8402974857881133471007083745436809127296054293775383549824742623937028497898215256929178577083970960121625602506027316549718402106494049978375604247408

24 12371712231207064758338744862673570832373041989012943539678727080484951695515930485641394550792153037191858028212512280926600304581386791094

23 5543542935523747700991431848906143793069037997096433133255695864648400840733488554456638692402087571124206008540851348293394572022 755970286667345339661519123315222619353103732072409481167391410479517925792743631234987038883317634987271171404439792

21 31374751050137102720420538137382214513103312193698723653061351991346433379389385793965576992246021316463868

20 3962892199823037560207299517133362502106339705739463771515237113377010682364035706704472064940398

19 1523344971704879993080742810319229690899454255323294555776029866737355060592877569255844

18 1782112840842065129893384946652325275167838065704767655931452474605826692782532

17 6344814611237963971310297540795524400449443986866480693646369387855336

16 68745445609149931587631563132489232824587945968099457285419306

15 2266745568862672746374567396713098934866324885408319028

14 227449714676812739631826459327989863387613323440

13 69450664761521361664274701548907358996488

12 64528039343270018963357185158482118

11 182413291514248049241470885236

10 1568758030464750013214100

9 41044208702632496804

8 3266598486981642

7 789360053252

6 575780564

5 1262816

4 8512

3 184

2 12

1201

おねえさん問題の答え(𝑛𝑛 = 1～26)

3

http://oeis.org/A007764/b007764.txt

4余白

5
乱択近似のアイデア
本研究:
マルコフ連鎖モンテカルロ(MCMC)法を用いた乱択近似数え上げアルゴリズムの提案．
記法の準備 Ω = {𝑠𝑠𝑠𝑠 単純経路}
 Ξ𝑘𝑘 = 長さ𝑘𝑘以下の𝑠𝑠𝑠𝑠 単純経路
𝑘𝑘 = 2𝑛𝑛, 2𝑛𝑛 + 2,2𝑛𝑛 + 4, … , 𝐿𝐿∗

ただし，𝐿𝐿∗ は最長路の長さ:
𝑛𝑛 + 1 2 − 1 (𝑛𝑛が奇数)∗𝐿𝐿 = �𝑛𝑛 + 1 2 − 2 (𝑛𝑛が偶数)
つまり，求めたいものは Ω (= |Ξ𝐿𝐿∗ |)．

自己帰着(self-reducibility)
Ξ𝐿𝐿∗
Ξ𝐿𝐿∗−2Ξ2𝑛𝑛+4 Ξ2𝑛𝑛+2Ξ𝐿𝐿∗ =⋅⋯⋅⋅ Ξ2𝑛𝑛Ξ𝐿𝐿∗−2 Ξ𝐿𝐿∗−4
Ξ2𝑛𝑛+2Ξ2𝑛𝑛
逆数を考えると...
アイデア 2.
Ξ𝑘𝑘
Ξ𝑘𝑘+2

アイデア 1.
Ξ2𝑛𝑛 = 最短経路の個数=2𝑛𝑛𝑛𝑛= [Ξ𝑘𝑘+2 の中で長さ𝑘𝑘以下の𝑠𝑠-𝑡𝑡単純経路の割合]

6
乱択近似アルゴリズム
パラメータ：
𝜏𝜏 （マルコフ連鎖の遷移回数），𝑀𝑀 （モンテカルロ法のサンプル）．
入力：𝑛𝑛 （格子のサイズ）．
出力：𝑍𝑍 （経路総数の近似値）．
Set 𝑍𝑍 ≔ 1;For (𝑘𝑘 = 2𝑛𝑛 + 2; 𝑘𝑘 < 𝐿𝐿∗ ; 𝑘𝑘 ≔ 𝑘𝑘 + 2){

Set 𝑋𝑋 ∈ Ξ𝑘𝑘 ; (𝑋𝑋はマルコフ連鎖MCの初期状態)
Set 𝑆𝑆 ≔ 0; (Sはカウンター)
for(i=0; i<M; i++){for(j=0; j<𝜏𝜏; j++){𝑋𝑋を更新(マルコフ連鎖の1回遷移)}}
if(𝑋𝑋 ∈ Ξ𝑘𝑘−2 ) S++;

Set 𝑍𝑍 ≔ 𝑍𝑍 ∗}Output 𝑍𝑍;

𝑀𝑀;𝑆𝑆

Ξ𝑘𝑘 からの一様サンプリングをしている

7
近似精度保証定理
任意の 𝜖𝜖 (0 < 𝜖𝜖 < 1)と任意の𝛿𝛿 (0 < 𝛿𝛿 < 1)に対して，
Ξ𝑘𝑘 (𝑘𝑘 = 2𝑛𝑛, 2𝑛𝑛 + 2, … , 𝐿𝐿∗ )からの一様標本の数を𝑀𝑀 = 12𝑛𝑛3 2𝑛𝑛2 𝜖𝜖 −1 2 ln 𝑛𝑛2 𝛿𝛿 −1 とすると，近似解𝑍𝑍はPr 1 − 𝜖𝜖 |Ω| ≤ 𝑍𝑍 ≤ 1 + 𝜖𝜖 |Ω| ≥ 1 − 𝛿𝛿を満たす．

8余白9マルコフ連鎖モンテカルロ法
目標Ξ𝑘𝑘 からの一様サンプリング

マルコフ連鎖モンテカルロ法のアイデア
1. 目的の分布を定常分布にもつマルコフ連鎖を設計．
2. 十分な回数推移させ、極限分布からサンプリング．

状態空間xxx・・・

xx

漸近的に定常分布

10

準備: s-t単純経路の表現s𝑠𝑠-𝑡𝑡単純経路ts
単連結な色塗りt

11
マルコフ連鎖の状態遷移
sscctt
マルコフ連鎖𝑀𝑀𝑀𝑀
状態空間: Ξ𝑘𝑘 ，

状態遷移: 𝑋𝑋 → 𝑋𝑋𝑋は以下のように定められる;
1. セルcを一様ランダムに選ぶ．

2. 状態𝑌𝑌を𝑋𝑋のセルcを反転したものとする．
3. 𝑌𝑌 ∈ Ξ𝑘𝑘 なら𝑋𝑋 ′ = 𝑌𝑌，そうでなければ𝑋𝑋 ′ = 𝑋𝑋とする．
Check!
• 𝑌𝑌(の境界線)は𝑠𝑠-𝑡𝑡経路• 長さ𝑘𝑘以下
実装上は高速化のため，もう少しべたに場合分けする．

12
極限分布
定理
マルコフ連鎖MCの極限分布はΞ𝑘𝑘 上の一様分布
証明の概略• マルコフ連鎖MCは既約 (Ξ𝑘𝑘 上で(強)連結)• マルコフ連鎖MCは非周期的• マルコフ連鎖MCは詳細釣合の式∀𝑋𝑋, 𝑌𝑌 ∈ Ξ𝑘𝑘 , Pr 𝑋𝑋 → 𝑌𝑌 = Pr(𝑌𝑌 → 𝑋𝑋)を満たす．

MCMC法設計の常套手段• 既約で非周期的な有限状態マルコフ連鎖は一意の極限分布をもつ．
• 詳細釣合の式∀𝑋𝑋, 𝑌𝑌 ∈ 𝑆𝑆, Pr 𝑋𝑋 → 𝑌𝑌 = Pr 𝑌𝑌 → 𝑋𝑋

を満たせば，定常分布は𝑆𝑆上の一様分布

マルコフ連鎖MCのmixing time(𝜏𝜏)は未解決

13
計算機実験
マルコフ連鎖の遷移回数
サンプリング回数𝜏𝜏 = 30

𝑀𝑀 = 107(近似解は5回の平均)

𝒏𝒏25

真の解*

近似解実行時間

8.40 * 10150

8.55 * 10150

1h 27munknown

2.09 * 𝟏𝟏𝟏𝟏𝟐𝟐𝟐𝟐𝟕𝟕

2h 4m

30

1.74 * 10163

50

unknown

100

unknown

200

unknown

26

1.78 * 10163

1h 33m6.35 * 𝟏𝟏𝟏𝟏𝟔𝟔𝟔𝟔𝟔𝟔

5h 44m6.07 * 𝟏𝟏𝟏𝟏𝟐𝟐𝟐𝟐𝟐𝟐𝟐𝟐

1.196 * 𝟏𝟏𝟏𝟏9667

23h 20m96h

実行時間(理論値)は O 𝑀𝑀𝑀𝑀𝑛𝑛2

14
おねえさんに教えてあげたい...
220.2415𝑛𝑛𝑛𝑛
経路の総数はおよそ10≃ 1.744 個(予想)

n と log10 𝑍𝑍の関係

格子数 n×n における経路数 𝑍𝑍：n = 1~26 の正しい経路数：n = 1~200の近似値：放物線 𝑦𝑦 = 0.2415𝑥𝑥 2

3予想𝑛𝑛2Ω ≥ 3log10 Ω𝑛𝑛2
120.301030 0.301030
2121.079181 0.269795
31842.264818 0.251646
48512
3.930032 0.245627
51262816
6.101340 0.244054
6 575780564
8.760257 0.243340
7 7.8936E+11 11.897275 0.242802
8 3.2666E+15 15.514096 0.242408
9 4.10442E+19 19.613252 0.242139
10 1.56876E+24 24.195556 0.241956
11 1.82413E+29 29.261056 0.241827
12 6.4528E+34 34.809748 0.241734
13 6.94507E+40 40.841676 0.241667
14 2.2745E+47 47.356885 0.241617
15 2.26675E+54 54.355403 0.241580
16 6.87454E+61 61.837244 0.241552
17 6.34481E+69 69.802419 0.241531
18 1.78211E+78 78.250935 0.241515
19 1.52334E+87 87.182798 0.241504
20 3.96289E+96 96.598012 0.241495
21 3.1375E+106 106.496580 0.241489
22 7.5597E+116 116.878505 0.241485
23 5.5435E+127 127.743787 0.241482
24 1.2372E+139 139.092430 0.241480
25 8.403E+150 150.924433 0.241479
Ωlog10 Ω15

3予想160.000000
140.000000

16𝑛𝑛2Ω ≥ 3

log10 Ωlog_10のプロット120.000000
100.000000
80.000000

log_1060.000000
40.000000
20.000000
0.000000

0.350000

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

log10 Ωのプロット
𝑛𝑛2 log_10/n^20.300000
0.250000
0.200000
log_10/n^20.150000
0.100000
0.050000
0.000000

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