情報基礎実習I （第１４回） 木曜４・５限 担当：北川 晃

プログラミング演習：複数のサイコロの出目 𝑛 𝑑 (≥2)個のサイコロを振ったときの， 出目の和を求めるプログラムを作れ． 𝑛 𝑑 =5 𝜆=12869 𝑐=6925 𝜇= 2 15 混合型合同法により，一様乱数を発生する． 𝑥 𝑛+1 = 𝜆 𝑥 𝑛 +𝑐 mod 𝜇 得られた乱数列を𝜇で割って，規格化乱数[0,1)を生成する． 規格化乱数に 7−1 =6をかけて1をかければ， 1≤𝑢<7なる一様乱数が求まる． 上記の乱数列の小数部分を切り捨てれば， サイコロの出目が得られる． これを 𝑛 𝑑 回繰り返し，和をとる． 試行を20回繰り返す

複数のサイコロの出目：プログラム例 Dim mu As Integer = 2 ^ 15 Sub Main() Dim x, n, d, s, nd As Integer n = 20 nd = 5 x = 1 For j As Integer = 1 To n s = 0 For i As Integer = 1 To nd x = uni_rnd(x) d = dice(x) s = s + d Next i Next j Console.WriteLine(s) End Sub 試行を𝑛回繰り返す サイコロ 𝑛 𝑑 個の目の和

複数のサイコロの出目：プログラム例 Function dice(x) Dim d As Integer d = Fix(6 * x / mu + 1) Return d End Function Function uni_rnd(x) Dim lambda, c As Integer lambda = 12869 c = 6925 x = (lambda * x + c) Mod mu Return x 初期値𝑥=1のとき， 21 19 14 22 12 18 9 20 16 8 … サイコロの目に変換 一様乱数のサブルーチン

横軸に0.1×𝑖× 𝑛 𝑑 ×6，縦軸に 𝑐 𝑖 をとって，CSV形式で書き出す． 複数のサイコロの出目の頻度 𝑛 𝑑 (≥2)個のサイコロの出目の和𝑠の頻度を調べよ． 出目の和𝑠が次のどの領域にあるか調べる． 𝑠≤ 𝑛 𝑑 ×6×0.05 𝑛 𝑑 ×6×0.05≤𝑠< 𝑛 𝑑 ×6×0.15 𝑛 𝑑 ×6×0.15≤𝑠< 𝑛 𝑑 ×6×0.25 𝑛 𝑑 ×6×0.25≤𝑠< 𝑛 𝑑 ×6×0.35 ⋮ 𝑛 𝑑 ×6×0.85≤𝑠< 𝑛 𝑑 ×6×0.95 𝑛 𝑑 ×6×0.95≤𝑠 カウンターを準備 𝑐 0 𝑐 1 𝑐 2 𝑐 3 ⋮ 𝑐 9 𝑐 10 横軸に0.1×𝑖× 𝑛 𝑑 ×6，縦軸に 𝑐 𝑖 をとって，CSV形式で書き出す．

複数のサイコロの出目の頻度：プログラム例 Dim mu As Integer = 2 ^ 15 Sub Main() Dim x, n, d, s, nd As Integer Dim Writer As New IO.StreamWriter( _ "c:\…\distribution_dice10.csv") Dim c(10) As Integer n = 100000 : nd = 5 x = 1 For i As Integer = 0 To 10 c(i) = 0 Next i 𝑛=100,000回程度 試行を繰り返す カウンターの初期値

複数のサイコロの出目の頻度：プログラム例 For j As Integer = 1 To n s = 0 For i As Integer = 1 To nd x = uni_rnd(x) d = dice(x) s = s + d Next i If s <= 6 * nd * 0.05 Then c(0) = c(0) + 1 ElseIf 6 * nd * 0.05 <= s And s < 6 * nd * 0.15 Then c(1) = c(1) + 1 ElseIf 6 * nd * 0.15 < s And s <= 6 * nd * 0.25 Then c(2) = c(2) + 1 ……… ElseIf 6 * nd * 0.95 < s Then c(10) = c(10) + 1 End If Next j 出た目に応じてカウンターの値を上げる

複数のサイコロの出目の頻度：プログラム例 For i As Integer = 0 To 10 'Console.WriteLine("c({0})={1}", i, c(i)) Writer.WriteLine("{0}, {1}", i, c(i)) Next i Writer.Close() End Sub 先のプログラミングのサイコロ乱数， 一様乱数の関数プロシージャを続ける

複数のサイコロの出目の頻度：出力結果 出目 平均値は 30+5 2 =17.5 最大値30 最小値 5

複数のサイコロの出目の平均，標準偏差 𝑛 𝑑 (≥2)個のサイコロを振ったときの， 出目の和の平均値𝐸と標準偏差𝜎を求めよ． 試行𝑗回目の出目の和を 𝑠 𝑗 とする（配列を用いる）． 平均値（期待値）𝐸を求める． 標準偏差𝜎を求める． 平均値𝐸，標準偏差𝜎の正規分布をプロットしてみよ． 𝑓 𝑥 = 1 2𝜋 𝜎 exp − 𝑥−𝐸 2 2 𝜎 2 𝐸= 1 𝑛 𝑗=1 𝑛 𝑠 𝑗 , 𝜎= 1 𝑛 𝑗=1 𝑛 𝑠 𝑗 −𝐸 2 0≤𝑥≤6 𝑛 𝑑

複数のサイコロの出目の頻度 Dim mu As Integer = 2 ^ 15 Sub Main() Dim x, n, d, s, nd As Integer Dim Writer As New IO.StreamWriter( _ "c:\…\distribution_dice10.csv") Dim Writer2 As New IO.StreamWriter( _ "c:\…\Gaussian.csv") Dim c(10) As Integer Dim d_sum(100000) As Integer, ave, sigma, x0, y0 As Single n = 100000 : nd = 5 x = 1 For i As Integer = 0 To 10 c(i) = 0 Next ガウス曲線の書き出し

各試行のさいころの出目の和を配列に書きだす 複数のサイコロの出目の頻度 For j As Integer = 1 To n s = 0 For i As Integer = 1 To nd x = uni_rnd(x) d = dice(x) 'Console.WriteLine(d) s = s + d Next i 'Console.WriteLine(s) If s <= 6 * nd * 0.05 Then c(0) = c(0) + 1 ……… ElseIf 6 * nd * 0.95 < s Then c(10) = c(10) + 1 End If d_sum(j) = s Next j 各試行のさいころの出目の和を配列に書きだす

複数のサイコロの出目の頻度 For i As Integer = 0 To 10 'Console.WriteLine("c({0})={1}", i, c(i)) Writer.WriteLine("{0}, {1}", i * 0.1 * nd * 6, c(i)) Next i Writer.Close() ave = ave_func(d_sum, n) Console.WriteLine(ave) sigma = sd_func(d_sum, ave, n) Console.WriteLine(sigma) For x0 = 0 To 6 * nd Step 0.1 y0 = 1 / Math.Sqrt(2 * Math.PI) _ * Math.Exp(-(x0 – ave) ^ 2 / 2 / sigma ^ 2) Writer2.WriteLine(“{0,12},{1,12}”, x0, y0) Next x0 Writer2.Close() End Sub 平均値の計算，書き出し 標準偏差の計算，書き出し ガウス曲線の書き出し

複数のサイコロの出目の頻度 Function ave_func(a, n) Dim s, ave As Single s = 0.0 For i As Integer = 1 To n s = s + a(i) Next i ave = s / n Return ave End Function Function sd_func(a, ave, n) Dim s, sd As Single s = s + (a(i) - ave) ^ 2 Next sd = Math.Sqrt(s / n) Return sd さらにサイコロ乱数， 一様乱数の関数を続ける．

複数のサイコロの出目の頻度：出力結果 3.8189 𝜎 𝐸 17.501 中心極限定理： 試行回数𝑛が十分大きい とき，正規分布に近づく 𝑓 𝑥 = 1 2𝜋 𝜎 exp − 𝑥−𝐸 2 2 𝜎 2 3.8189 𝜎 𝐸 17.501

プログラミング演習:正規乱数の発生 出現確率が正規分布に従うような乱数列 （正規乱数）を生成するプログラムを作れ． 正規一様乱数[0,1)の𝑘個の平均をとる 𝐸= 1 𝑛 𝑖=1 𝑛 𝑢 1 + 𝑢 2 +⋯+ 𝑢 𝑘 𝑘 = 1 2 平均値： 𝜎 2 = 1 𝑛 𝑖=1 𝑛 𝑥 𝑖 −𝐸 2 → 0 1 𝑥− 1 2 2 𝑑𝑥 = 1 12 分散： 𝑛→大のとき，平均値𝐸，分散 𝜎 2 𝑛 の正規分布 𝑧 𝑖 = 𝑥 𝑖 − 1 2 1 12𝑛 乱数　　　　　　　は，規格化正規分布に近づく

プログラミング演習:正規乱数の発生 𝐸= 1 𝑛 𝑖=1 𝑛 𝑢 1 + 𝑢 2 +⋯+ 𝑢 𝑘 𝑘 = 1 2 正規一様乱数[0,1)の𝑘個の平均をとる 𝐸= 1 𝑛 𝑖=1 𝑛 𝑢 1 + 𝑢 2 +⋯+ 𝑢 𝑘 𝑘 = 1 2 平均値： 𝜎 2 = 1 𝑛 𝑖=1 𝑛 𝑥 𝑖 −𝐸 2 → 0 1 𝑥− 1 2 2 𝑑𝑥 = 1 12 分散： 𝑘→大のとき，平均値𝐸，分散 𝜎 2 𝑘 の正規分布 𝑧 𝑖 = 𝑥 𝑖 − 1 2 1 12𝑘 乱数　　　　　　　は，規格化正規分布に従う 𝑓 𝑧 = 1 2𝜋 exp − 𝑧 2 2 −5.5≤𝑧<5.5を，幅Δ𝑧=1の区間に分割し，頻度を調べよ．

正規乱数の発生：出力例 𝑓 𝑧 = 1 2𝜋 exp − 𝑧 2 2

正規乱数の発生：プログラム例 Dim mu As Integer = 2 ^ 15 Sub Main() Dim k, n, x, c(10) As Integer Dim z, s, p(10) As Single Dim Writer As New IO.StreamWriter( _ "d:\…\random_number_normal.csv") k = 5 : n = 10000 x = 1 For i As Integer = 1 To 10 c(i) = 0 Next i

正規乱数の発生：プログラム例 For i As Integer = 1 To n s = 0.0 For j As Integer = 1 To k x = uni_rnd(x) s = s + x / mu Next j s = s / k z = (s - 1 / 2) / Math.Sqrt(1 / 12 / k) Console.WriteLine(z) If z <= -4.5 Then c(0) = c(0) + 1 ElseIf -4.5 < z And z <= -3.5 Then c(1) = c(1) + 1 ……… ElseIf 4.5 < z And z <= 5.5 Then c(10) = c(10) + 1 End If Next i 平均値の計算 変数の変換

正規乱数の発生：プログラム例 For i As Integer = 0 To 10 p(i) = c(i) / n Writer.WriteLine("{0}, {1}", i - 5, p(i)) Next Writer.Close() End Sub Function uni_rnd(x) Dim lambda, c As Integer lambda = 12869 c = 6925 x = (lambda * x + c) Mod mu Return x End Function 変数𝑧の頻度の書き出し

規格化正規分布曲線のためのプログラム Sub Main() Dim Writer As New IO.StreamWriter("d:\…\Gaussian_curve.csv") Dim y, E, sigma As Single E = 0 : sigma = 1 For x As Single = -5 To 5 Step 0.1 y = Gauss(x, E, sigma) Writer.WriteLine("{0,12}, {1,12}", x, y) Next x Writer.Close() End Sub Function Gauss(x, E, sigma) Dim y As Single y = 1 / Math.Sqrt(2 * Math.PI * sigma ^ 2) _ * Math.Exp(-(x - E) ^ 2 / 2 / sigma ^ 2) Return y End Function

プログラミング演習：指数乱数列 乱数列の発生頻度が に従うものを，指数乱数列という． 𝑓 𝑦 =𝑎 𝑒 −𝑎𝑦 　乱数列の発生頻度が　　　　　　　　　 に従うものを，指数乱数列という． 規格化一様乱数列 𝑢 1 , 𝑢 2 ,⋯, 𝑢 𝑛 から指数乱数列 𝑦 1 , 𝑦 2 ,⋯, 𝑦 𝑛 を 生成するプログラムを作れ． 𝑓 𝑦 =𝑎 𝑒 −𝑎𝑦 𝑦≥0, 𝑎>0 逆関数法 𝐹 𝑦 = 0 𝑦 𝑎 𝑒 −𝑎𝑥 𝑑𝑥 =𝑎 1− 𝑒 −𝑎𝑦 0≤𝐹 𝑦 ≤1であり，規格化一様乱数 𝑢 𝑖 より 𝑦 𝑖 =− 1 𝑎 log 1− 𝑢 𝑖 は 指数乱数を与える． 出現頻度を0≤𝑦≤10の領域について計算し，分布関数を重ねてみよ．

指数乱数の発生：出力例 𝑓 𝑦 = 𝑒 −𝑦

指数乱数の発生：プログラム例 Dim mu As Integer = 2 ^ 15 Sub Main() Dim Writer As New IO.StreamWriter( _ "d:\…\randan_number_exponential.csv") Dim x, n, c(10) As Integer Dim a, y, p(10) As Single n = 10000 a = 1 x = 1 For i As Integer = 1 To 10 c(i) = 0 Next

指数乱数の発生：プログラム例 For i As Integer = 1 To n x = uni_rnd(x) y = -1 / a * Math.Log(1 - x / mu) If y >= 0 And y < 1.0 Then c(1) = c(1) + 1 ElseIf y >= 1.0 And y < 2.0 Then c(2) = c(2) + 1 ……… ElseIf y >= 8.0 And y < 9.0 Then c(9) = c(9) + 1 Else c(10) = c(10) + 1 End If Next

指数乱数の発生：プログラム例 For i As Integer = 1 To 10 p(i) = c(i) / n Writer.WriteLine("{0}, {1}", i - 0.5, p(i)) Next Writer.Close() End Sub Function uni_rnd(x) Dim lambda, c As Integer lambda = 12869 c = 6925 x = (lambda * x + c) Mod mu Return x End Function