有一种赌博游戏,投掷一枚硬币,如果正面向上,赌徒将获得参赌赌资的1+a倍,如果反面向上,赌徒将获得参赌赌资的1-b倍,赌局连续进行。现在利用蒙特卡洛模拟方法来说明ab取值对赌局结果的影响。
为方便起见,设置赌徒的初始资金为1 a=0.2。
a=0.2
假设连续进行M次赌博,b取值为0.2,M取值为1000。
def dubo(M,b):
rec=np.ones(M)
for dujv in range(1,M):
if np.random.rand()>0.5:
rec[dujv]=rec[dujv-1]*(1+a)
else:
rec[dujv]=rec[dujv-1]*(1-b)
return rec
res=dubo(1000,0.2)
def showline(sz):
plt.figure(figsize=(16,9),frameon=False,dpi=120)
plt.plot(sz,'.-')
plt.xlabel('赌局次数')
plt.ylabel('赌徒持有赌资')
plt.show()
showline(res)
通过这个赌资随赌局变化图上看出,虽然赌资最大值达到8,赌徒进行1000次赌博之后赌资趋近于0,在大约前400次赌博之后即可认为该赌徒已经彻底输光了。
如果减小b的值,看一下情况会不会有所变化。
取b=0.15
res2=dubo(1000,0.15)
showline(res2)
当b取0.15时,赌徒的赌资迅速增加,等到进行后面几轮时,赌资已经达到了惊人的$10^{6}$数量级,这时候可以认为赌场已经实质性破产了。
从b=0.2时赌徒输光到b=0.15赌场破产,在其中会有一个过渡。现在开始寻找这个临界值。
由于每次赌博的赌资变化有着不确定性,为了能够反映当赌局次数趋向于正无穷的情况,可以进行多次模拟并取平均。为了简化问题,只对每次模拟最终赌资是否大于初始赌资,如果模拟结束时最终赌资大于初始赌资,则赌徒获胜。通过统计不同b取值的情况下赌徒的胜率,找出b的临界点。
为研究问题方便,每个b的取值进行5000次模拟,每次模拟进行1000轮赌局。
Mt=5000
Lt=1000
def dubot(b):
r=np.zeros_like(b)
for bdog in range(0,len(b)):
ui=np.zeros(1)
bb=b[bdog]
for cdog in range(0,Mt):
uii=np.ones(1)
for ldog in range(0,Lt):
if np.random.rand()>0.5:
uii=uii*(1+a)
else:
uii=uii*(1-bb)
if (uii>1):
ui=ui+1
r[bdog]=ui/Mt
return r
bbb2=np.linspace(0.1,0.2,21)
rrr2=dubot(bbb2)
可以很轻易的画出赌徒胜率与b的关系图
plt.figure(figsize=(15,9),frameon=False,dpi=120)
plt.plot(bbb2,rrr2,'x-')
plt.xlabel('b')
plt.ylabel('赌徒胜率')
# plt.yticks(np.arange(0,1,0.1))
# plt.xticks(np.arange(0.15,0.2,.005))
plt.show()
从图中可以看出来大约在0.16到0.18之间存在着一个胜率为0.5的临界点,可以继续放大图像查看。
bbb=np.linspace(0.15,0.2,101)
rrr=dubot(bbb)
plt.figure(figsize=(15,9),frameon=False,dpi=120)
plt.plot(bbb,rrr,'.-')
plt.xlabel('b')
plt.ylabel('赌徒胜率')
plt.yticks(np.arange(0,1,0.1))
plt.xticks(np.arange(0.15,0.2,.005))
plt.show()
从图中可以看出随着b的增大,赌徒的胜率逐渐减少以至于趋向于0,当b接近0.167时胜率接近于0.5。
从理论上来说,b=0.167作为临界相变点可以从细致平衡上面来理解。此时b的值满足关系$(1-b)(1+a)=1$即在第一次输钱之后,第二次再赢钱后的结果将回到第一次输钱之前。如果b值大于临界值,在经过足够多的赌博轮数之后可以认为是必输的,至少游戏结束时比初始钱数大的可能性几乎为0。
当b值大于临界值的情况下,寻求进行赌局轮数的最优解。
现在取b=0.19的情况,进行五万次模拟,每次模拟2000轮。
在b等于0.19的情况下,如果仅仅考虑数学期望,这是大于1的,按道理应该是应该一直玩下去的,但实际上虽然期望值大于1
但达到这个结果却并不容易
MM=50000
NN=2000
quu=np.ones([MM,NN])
for ui in range(1,NN):
quu[:,ui]= quu[:,ui-1]*(np.random.randint(2,size=MM)*.39+0.81)
首先看一下如果进行n次游戏后的赌资的平均值。
ssl2=np.mean(quu,axis=0)
plt.figure(figsize=(15,9),frameon=False,dpi=120)
plt.plot(ssl2,'.-')
plt.xlabel('次数')
plt.ylabel('赌资平均数')
plt.show()
通过对这50000次的赌博进行平均,可以发现平均值都是大于1,完全可以一直玩下去。
看一下中位数,马上就会发现问题,所有这五万次的模拟表明中位数全部小于1,即会有大于一半的可能性会输钱。
ssl=np.median(quu,axis=0)
plt.figure(figsize=(15,9),frameon=False,dpi=120)
plt.plot(ssl,'.-')
plt.xlabel('次数')
plt.ylabel('赌资中位数')
plt.show()
为了更直观的看出分布情况,现在截取第125次、250次、500次、1000次结束后的直方分布图。
plt.figure(figsize=(15,9),frameon=False,dpi=120)
sj=plt.hist(quu[:,125],log=True,bins=np.linspace(0,10,101))
plt.xlabel('赌资')
plt.ylabel('频数')
plt.show()
第125次
plt.figure(figsize=(15,9),frameon=False,dpi=120)
sj=plt.hist(quu[:,250],log=True,bins=np.linspace(0,10,101))
plt.xlabel('赌资')
plt.ylabel('频数')
plt.show()
第250次
plt.figure(figsize=(15,9),frameon=False,dpi=120)
sj=plt.hist(quu[:,500],log=True,bins=np.linspace(0,10,101))
plt.xlabel('赌资')
plt.ylabel('频数')
plt.show()
第500次
plt.figure(figsize=(15,9),frameon=False,dpi=120)
sj=plt.hist(quu[:,1000],log=True,bins=np.linspace(0,10,101))
plt.xlabel('赌资')
plt.ylabel('频数')
plt.show()
第1000次
从中可以明显看出来随着次数的增加,分布会逐渐向0靠拢,即全部赔光的可能越来越多。
到现在可以得出结论,尽管单次赌博的期望值是大于1的,但从统计上来看随着次数的增加,赢钱的可能性越来越小,虽然一旦赢钱所得到的数字会很大以至于将整体的期望值提高到1以上。归根结底原因在于,当b大于临界点的数值情况下,减少一次后需要增加一次以上才能回到初始值,而增加与减少是等概率的,就是这个原因造成了虽然期望值是大于1的,但实际操作起来钱却越来越少。
所以得出的结论是如果b大于临界点的数值,就不要进行连续游戏。如果期望值大于1,则可以进行单次的不会进行连续累计的游戏。