サイコロ振って中心極限定理を確かめる

 中心極限定理とは、

中心極限定理は標本平均と母平均との誤差を論ずるものである。多くの場合、母集団の分布がどんな分布であっても、その誤差は標本の大きさを大きくしたとき近似的に正規分布に従う。wikipedia

とあるように標本を大きくしていくと正規分布に近づくという定理です。

サイコロを想定してpythonで表現してみます。

numpyのrandintでランダムに1〜6の目を出すサイコロを回数分出します。

とりあえず10000回くらい振ってみましょう。

view plainprint?

1. #! /usr/bin/env python  
2.   
3. # -*- coding:utf-8 -*-  
4.   
5. import pandas as pd  
6. import numpy as np  
7.   
8.   
9. print("[dice]")  
10. print()  
11.   
12. # サイコロの目(1〜6)をランダムにnumpy配列で出力  
13. pip = np.random.randint(1, 7, 10000)  
14.   
15. # 作られた配列の先頭を表示  
16. print(pip[:10])  
17. print()  
18. print(f"平均値: %f" % pip.mean())  
19. print()  

#! /usr/bin/env python

# -*- coding:utf-8 -*-

import pandas as pd
import numpy as np


print("[dice]")
print()

# サイコロの目(1〜6)をランダムにnumpy配列で出力
pip = np.random.randint(1, 7, 10000)

# 作られた配列の先頭を表示
print(pip[:10])
print()
print(f"平均値: %f" % pip.mean())
print()

そしてmatplotlibのヒストグラムで分布を見てみましょう。

view plainprint?

1. #ヒストグラムでプロット  
2. import matplotlib.pyplot as plt  
3. import japanize_matplotlib  
4.   
5. plt.title("サイコロ")  
6. plt.xlabel("出た目")  
7. plt.ylabel("出た回数")  
8. plt.hist(pip, bins=6, ec="k")  
9. plt.show()  

#ヒストグラムでプロット
import matplotlib.pyplot as plt
import japanize_matplotlib

plt.title("サイコロ")
plt.xlabel("出た目")
plt.ylabel("出た回数")
plt.hist(pip, bins=6, ec="k")
plt.show()

平均値: 3.502500

1~6までほぼ同じ回数出てるのが分かります。

期待値は(1+2+3+4+5+6) / 6 = 7/2 = 3.5なので平均値もほぼ合ってます。

サイコロ1つなら感覚的にもつかみやすいでしょう。

次に2個のサイコロを振ってみます。

view plainprint?

1. import pandas as pd  
2. import numpy as np  
3. import matplotlib.pyplot as plt  
4. import japanize_matplotlib  
5. #最瀕値を出すために使用  
6. import statistics  
7.   
8. #振る回数  
9. roll = 10000  
10.   
11. #1つ目のサイコロ  
12. pip = np.random.randint(1, 7, roll)  
13. #2こ目のサイコロの  
14. pip2 = np.random.randint(1, 7, roll)  
15.   
16. print("[dice]")  
17. print()  
18. #先頭10回分の目  
19. print(pip[:10])  
20. print(pip2[:10])  
21.   
22. sm = pip + pip2  
23. print()  
24. print(f"平均値: %f" % sm.mean())  
25. print(f'最瀕値: %f' % statistics.mode(sm))  
26.   
27. plt.title("サイコロ")  
28. plt.xlabel("出た目")  
29. plt.ylabel("出た回数")  
30. plt.hist(sm, bins=11, ec="k")  
31.   
32. plt.show()  

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import japanize_matplotlib
#最瀕値を出すために使用
import statistics

#振る回数
roll = 10000

#1つ目のサイコロ
pip = np.random.randint(1, 7, roll)
#2こ目のサイコロの
pip2 = np.random.randint(1, 7, roll)

print("[dice]")
print()
#先頭10回分の目
print(pip[:10])
print(pip2[:10])

sm = pip + pip2
print()
print(f"平均値: %f" % sm.mean())
print(f'最瀕値: %f' % statistics.mode(sm))

plt.title("サイコロ")
plt.xlabel("出た目")
plt.ylabel("出た回数")
plt.hist(sm, bins=11, ec="k")

plt.show()

平均値: 7.044900

最瀕値: 7.000000

サイコロを2つにすると足した7が一番多い山のような形になることが分かります。

振るサイコロ6個の目の合計ではどうでしょう？

サイコロ6つだと、出る目は6の6乗の46656とおりになるので、それ以上回数を振らないとまばらな結果になりそうです。

10万回くらい振ってみましょう

view plainprint?

1. #! /usr/bin/env python  
2.   
3. # -*- coding:utf-8 -*-  
4.   
5. #  
6. import pandas as pd  
7. import numpy as np  
8. import matplotlib.pyplot as plt  
9. import japanize_matplotlib  
10. #最瀕値を出すために使用  
11. import statistics  
12.   
13. #振る回数  
14. roll = 100000  
15.   
16. # ndarrayの初期化  
17. sm = np.zeros(roll)  
18.   
19. print("[dice]")  
20.   
21.   
22. for i in range(6):  
23.     pip = np.random.randint(1, 7, roll)  
24.       
25.     #出ための先頭10回  
26.     print(pip[:10])  
27.     sm += pip  
28.   
29. print()  
30. print(f"平均値: %f" % sm.mean())  
31. print(f"最瀕値: %f" % statistics.mode(sm))  
32. print()  
33. #出ための合計の種類  
34. print(np.unique(sm))  
35.   
36.   
37. plt.figure(figsize=(10, 8))  
38. plt.title("サイコロ")  
39. plt.xlabel("出た目の合計")  
40. plt.ylabel("出た回数")  
41. plt.grid()  
42. plt.hist(sm, bins=31, ec="k")  
43.   
44. plt.show()  

#! /usr/bin/env python

# -*- coding:utf-8 -*-

#
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import japanize_matplotlib
#最瀕値を出すために使用
import statistics

#振る回数
roll = 100000

# ndarrayの初期化
sm = np.zeros(roll)

print("[dice]")


for i in range(6):
    pip = np.random.randint(1, 7, roll)
    
    #出ための先頭10回
    print(pip[:10])
    sm += pip

print()
print(f"平均値: %f" % sm.mean())
print(f"最瀕値: %f" % statistics.mode(sm))
print()
#出ための合計の種類
print(np.unique(sm))


plt.figure(figsize=(10, 8))
plt.title("サイコロ")
plt.xlabel("出た目の合計")
plt.ylabel("出た回数")
plt.grid()
plt.hist(sm, bins=31, ec="k")

plt.show()

平均値: 21.013420

最瀕値: 21.000000

ほぼ正規分布に近づきました。

標本を大きくすれば正規分布に近づいていくという中心極限定理が感覚的につかめました。