利用python和matplotlib对囚犯与盒子问题做可视化分析

几天前在看到一道有趣的谜题
: The condemned prisoners and the boxes
，作者给出了解法和并使用
R
来展示分析过程，使用
data visualization
来分析理解问题的方法很赞。

 

R
我没有用过，作为一个
python
爱好者，相信类似的统计分析和可视化
python
一定是有类似的工具的，正好几天前工作中需要分析并可视化数据，就开始寻找，最终在
stackoverflow
上发现很多人推荐
matplotlib
这个工具，配合
ipython
使用，体验很不错，所以遇到这道谜题的时候，用
python+matplotlib
也实验了一把，谜题背后的原理和数据证明参考最后给出的链接，这里主要通过例子来介绍如何使用
python
来做可视化。

 

进入正题：

有
100
个罪犯，编号
0-99
，
又有
100
个盒子，每个盒子里放一张带编号的纸条，也是
0-99
，对应
100
个罪犯，盒子是完全一样的。这
100
个盒子随机排列，成一条直线放在一个房间中。

游戏开始，第一个囚犯进入这个房间，选择
50
个盒子打开，如果有对应自己编号的纸条，视为成功，否则失败，然后离开房间到另外一个封闭的房间，盒子全部合上恢复原样。然后第二个囚犯进入，同样选择
50
个盒子打开，依次类推；

游戏开始之前囚犯可以商量制定策略，开始后囚犯没有机会再相互沟通；

如果
100
个囚犯都找到了自己的编号，这些囚犯就算赢了，可以释放。但是如果有任何一个囚犯没有找到自己的编号，则所有囚犯都不能释放；

囚犯应该采取什么策略才能使大家能被释放的几率最大？

 

简单考虑一下这个问题，能得到几点思路：

1.      

如果所有囚犯都随机挑选
50
个盒子，每个人找到自己编号的概率都为
0.5
，
100
个囚犯都找到自己盒子的概率就是
0.5**100,
几乎为
0
，所以随机选择可以说必输无疑
;

2.      

不管策略如何，由于盒子的排列完全随机，第一个囚犯不管使用什么策略，获胜的几率恒为
0.5
，所以所有囚犯都选中的概率不可能大于
0.5;

3.      

由于策略制定后才开始游戏，囚犯之间不能再交流信息，所以不管采用什么策略，对任何一个单个的囚犯来说，由于没有其他信息，箱子的排列还是完全随机的，那么他选中的概率也恒为
0.5
，也就是说，所有的囚犯，他能够选中的概率都是
0.5
，如果把这些随机事件的概率表示为
P0-P99,
则
P0==P1==…==P99==0.5
；

4.      

当
P0,…,P99
不相关时，猜中的囚犯的数目
PN(n)
服从二项分布，
PN(100)=0.5**100
，几乎为
0
。显然制定的策略需要使
P0,…,P99
是相关的；

5.      

一个好的策略，应该尽量使
P0-P99
“完全相关”，即要错大家绝大部分都猜错，要猜对大家全部都猜对，这样才能使
PN(100)
尽量接近上限
0.5

 

数学牛的同学当然可以用纯数学的方法分析解决这个问题，但是最简单直接的方法就是利用
python
程序来做随机模拟了
:)

 

首先需要定义一个模拟框架，来分析不同的策略会产生什么样的效果：

1.      

盒子的排列可以抽象成一个数组
boxes
，比如
[4,0,2,1,3]
表示
5
个盒子从左到右分别放了标签
4
，
0
，
2
，
1
，
3

2.      

一个策略以一个函数
strategy(boxes, p)
表示，
p
为囚犯的标签，函数表示标签为
p
的囚犯对
boxes
的选择策略，返回
p
是否成功找到了自己的标签。

 

框架代码如下：

import sys
import random
def strategy(boxes, p):
''' return True if strategy success '''
return False
def simulate(n, strategy, times):
''' n: number of boxes or prisoners
strategy: strategy used
times: random simulation times
return numbers of succeeded prisoners as list
'''
boxes = range(n)
result = []
for i in xrange(times):
random.shuffle(boxes)
success = 0
for p in xrange(n):
if strategy(boxes, p):
success += 1
result.append(success)
return result

 

先来定义一个随机策略，看看效果如何：

 

def random_strategy(boxes, p):
''' pick half of boxes at random '''
choices = range(len(boxes))
random.shuffle(choices)
picks = choices[:len(boxes)/2]
for i in picks:
if boxes[i] == p:
return True
return False
result = simulate(100, random_strategy, 100)
print result

 

结果为：

[51, 55, 50, 51, 53, 53, 54, 56, 58, 51, 59, 44, 45, 53, 54, 53, 49, 49, 51, 48, 52, 52, 51, 47, 54, 50, 53, 57, 50, 47, 42, 45, 49, 51, 48, 51, 49, 62, 56, 51, 51, 45, 48, 54, 54, 48, 53, 44, 59, 55, 43, 50, 57, 47, 44, 50, 51, 48, 47, 44, 53, 53, 53, 36, 52, 43, 45, 59, 53, 48, 59, 55, 50, 48, 43, 54, 45, 57, 50, 49, 45, 46, 56, 52, 55, 46, 56, 55, 41, 54, 58, 54, 47, 47, 51, 48, 40, 49, 50, 50]

 

从数据上看符合符合二项分布，近似符合正态分布，是不是这样呢，用
matplotlib
来可视化数据：

 

import matplotlib.pyplot as plt
plt.hist(result)
plt.show()

 



 

 

 

有点正太分布的雏形了，不过样本数据少，图形不太漂亮
...,

增加模拟次数再试试：

 
 

n = 100
result = simulate(n, random_strategy, 4000)
print result
import matplotlib.pyplot as plt
dist = [result.count(i) for i in range(n+2)]
plt.bar(range(n+2), dist)
plt.show()

 
 



 

 

 

漂亮多了吧
:)

 

随机策略就分析到这里，下面先来看看标准解法：

 

每个囚犯
p

按照如下序列选择盒子，直到选到盒子内的标签为
p

止：

p, boxes[p], boxes[boxes[p]], boxes[boxes[boxes[p]]] …

比如
boxes=[3,0,1,4,2],

选择顺序为：

p=0: 
3, 4, 2, 1, 0

p=1: 
0, 3, 4, 2, 1

 

simple & elegant

！

用代码表示：

 

 

def standard_strategy(boxes, p):
''' pick boxes start with p, then boxes[p], then boxes[boxes[p]] ...
until we find p '''
times_remain = len(boxes)/2
current = p
while times_remain > 0:
times_remain -= 1
if boxes[current] == p:
return True
else:
current = boxes[current]
return False
n = 100
result = simulate(n, standard_strategy, 4000)
print result
import matplotlib.pyplot as plt
dist = [result.count(i) for i in range(n+2)]
plt.bar(range(n+2), dist)
plt.show()

  



 

相比随机策略几乎没可能全部猜对的绝望，这种策略让所有囚犯能都猜对的概率提升到了
12xx/4000,

大概
1/3

！

 

概率分布图也很有意思，一点不像概率论上学到的经典概率分布，不是连续平滑分布的；

猜中的囚犯数为
0-50

有大概
2/3

的概率，为
51-99

的概率为
0

，为
100

的概率大概
1/3,

这貌似反映了“要错大家绝大部分都猜错，要猜对大家全部都猜对”的思路。

 

至于谜题背后的数学原理以及
matplotlib

安装，可以参考下面几个链接：

  http://matplotlib.sourceforge.net/users/installing.html
 matplot的安装

  http://www.statisticsblog.com/2010/07/100-prisoners-100-lines-of-code
 最开始提到的blog post

  http://www.mast.queensu.ca/~peter/inprocess/prisoners.pdf
 对这个问题的数学分析

  http://mathworld.wolfram.com/PermutationCycle.html
 PermutationCycle的相关资料