最終回　乱数思考

宝クジの購入は非合理的

人間は確率的な思考が得意ではありません。宝クジが当たる確率は非常に低く、賞金の期待値は投資額の半分しかありません。宝クジの購入が非合理的なのは明らかです。にもかかわらず、たくさんの人が一攫千金を夢見て購入するのは不思議です。

宝クジの場合は嫌なら買わなければよいのですが、保険の種類の選択に悩んだり、携帯電話の料金プランに悩んだりと、確率計算と無縁に生活することは困難です。

間違えやすい確率問題

確率計算が必要な問題は、直感に反することがよくあるので注意が必要です。

誕生日一致問題

確率の見積りを間違えやすい例として、「⁠N人の人間がいるとき、同じ誕生日の人がいる確率はどれぐらいか？」という問題があります。誕生日は365種ありますし、自分と同じ誕生日の人を知っていることは珍しいので、かなり多くの人間を集めないと誰かの誕生日が一致することはないだろうと思いがちです。しかし実際は、50人ぐらい人間がいれば、ほぼ確実に誕生日が同じである人がいるはずです。

数式による計算

誕生日が一致する2人がいる確率の計算は比較的簡単です。N人の人間の誕生日がすべて異なる確率を計算し、これを1から引けばいいわけです。全員を並べたとき最初の人と2番目の人の誕生日が異なる確率は364/365であり、3番目の人の誕生日がまた異なる確率は363/365であり……ということになるので、N人の誕生日がすべて異なる確率は、

p = 1 * 364/365 * 363/365 * ... (365+1-N)/365

という式で計算できます。

50人いた場合に誕生日が一致する2人がいる確率（1-p）は、次のRubyプログラムで計算できます。

n = (ARGV.shift || 50).to_i

p = 1.0
n.times { |i|
  p *= (365.0 - i) / 365.0
}
puts 1.0 - p

実行結果は次のようになります。

% ruby tanjoubi0.rb 50
0.97037357

50人いた場合に誕生日が一致する2人がいる確率は、約97％となりました。この結果は正しいのでしょうか。

シミュレーションによる計算

数式による計算結果が心配な場合、シミュレーションによる計算も行えば、結果を確信できるでしょう。N人の誕生日をランダムに生成し、同じ誕生日の人が見つかるかどうかを何度も計算します。

n = (ARGV.shift || 50).to_i # 人間の数
trials = (ARGV.shift || 10000).to_i # 試行回数

match = 0 # 同じ誕生日の人がいた回数

trials.times {
  birthday = {}
  n.times { |i|
    day = rand(365)
    if birthday[day] then # 誕生日の衝突発見
      puts "#{birthday[day]}番目と#{i}番目が一致"
      match += 1
      break
    end
    birthday[day] = i
  }
}

puts match.to_f / trials.to_f

実行結果は次のようになります。

% ruby tanjoubi.rb 50
0番目と20番目が一致
2番目と22番目が一致
...
15番目と23番目が一致
5番目と26番目が一致
0.9728

50人いた場合に誕生日が一致する2人がいる確率も、約97％となりました。シミュレーションなので毎回結果は異なりますが、何度計算しても同じような値が得られるので、数式による計算が正しかったことをほぼ確信できるでしょう。

兄弟の性別問題

誕生日一致問題の場合、きちんとした計算が必要であることはわかりますが、直感が強力な場合は間違った答えをなかなか訂正できないことがあります。

「2人の子どもの一方が男の子の場合、もう一方も男の子である確率は？」という問題も間違えやすいことで有名です。この確率は1/2だと考えてしまいがちですが、正解は1/3です。

場合分けによる計算

2人の子どもをA、Bとするとき、A、Bが男の子か女の子かという場合分けは図1のようになります。

一方が男の子であるケースは3種類ありますが、そのうち2種類のケースでもう一方は女の子になるので、この問題の答えは1/3だということになります。

シミュレーションによる計算

この問題の場合も、実際にシミュレーションによる計算をしてみることで答えを確信できます。

trials = (ARGV.shift || 10000).to_i # 試行回数
boys = 0; girls = 0
BOY = 0; GIRL = 1

trials.times {
  # 2人の子どもの性別をランダムにセット
  gender = []
  gender[0] = rand(2)
  gender[1] = rand(2)

  # 一方が男の子だった場合だけ計算
  if gender[0] == BOY || gender[1] == BOY then
    # 男の子のどちらかを選択してind1に格納
    if gender[0] == BOY && gender[1] == BOY
      ind1 = rand(2)
    elsif gender[0] == BOY
      ind1 = 0
    elsif gender[1] == BOY
      ind1 = 1
    end
    # もう一方をind2に格納
    ind2 = (ind1 == 0 ? 1 : 0)

    boys += 1 if gender[ind2] == BOY
    girls += 1 if gender[ind2] == GIRL
  end
}
puts "BOY: #{boys}"
puts "GIRL: #{girls}"
puts "Ratio: #{boys.to_f/(boys+girls)}"

実際に男の子／女の子をランダムに計算し、どちらかが男の子だったときにもう一方の性別を調べます。2人とも男の子だったとき、どちらを選んだ場合でももう一方が男の子であることは確かなのですが、問題の条件に忠実に従うように、ランダムに一方を選んでからもう一方の性別を調べるようにしています。

実行結果は次のようになります。

% ruby danshi.rb
BOY: 2487
GIRL: 5023
Ratio: 0.33115845539281

シミュレーションでも約33％（1/3）という答えを得ることができました。

火曜日に生まれた男の子問題

先の例を少し複雑にしたものとして、「⁠2人の子どもの一方が男の子で火曜日生まれの場合、もう一方も男の子である確率は？」という問題を考えてみましょう。「⁠火曜日生まれ」という条件が意味不明なので、この問題の答えは1/2か1/3だと考えるのが普通ですが、あらゆる可能性をきちんと考えて計算すると確率は13/27=0.48148になることがわかります。

実際にシミュレーションするRubyプログラムは次のようになります。

trials = (ARGV.shift || 100000).to_i # 試行回数
boys = 0; girls = 0
BOY = 0; GIRL = 1

trials.times {
  # 子どもの性別と生まれた曜日をランダムにセット
  gender = []
  gender[0] = rand(2)
  gender[1] = rand(2)
  weekday = []
  weekday[0] = rand(7)
  weekday[1] = rand(7)

  # 一方が男の子で火曜日生まれだった場合だけ計算
  if (gender[0] == BOY && weekday[0] == 2) ||
      (gender[1] == BOY && weekday[1] == 2) then
    # 男の子のどちらかを選択してind1に格納
    if (gender[0] == BOY && weekday[0] == 2) &&
        (gender[1] == BOY && weekday[1] == 2) then
      ind1 = rand(2)
    elsif gender[0] == BOY && weekday[0] == 2 then
      ind1 = 0
    elsif gender[1] == BOY && weekday[1] == 2 then
      ind1 = 1
    end
    # もう一方をind2に格納
    ind2 = (ind1 == 0 ? 1 : 0)
    boys += 1 if gender[ind2] == BOY
    girls += 1 if gender[ind2] == GIRL
  end
}

puts "BOY: #{boys}"
puts "GIRL: #{girls}"
puts "Ratio: #{boys.to_f/(boys+girls)}"

実行結果は次のようになります。

% ruby kayou.rb
BOY: 6619
GIRL: 7233
Ratio: 0.477837135431707

確率は1/2より少しだけ小さな値になることがわかります。

モンティ・ホール問題

確率の見積りを誰もが間違えやすい問題として、「⁠モンティ・ホール問題」が有名です。これはモンティ・ホールという司会者のゲームショー番組に由来するもので、次のようなものです。

「プレイヤーの前に3つのドアがあって、1つのドアの後ろには景品の新車が、2つのドアの後ろにはヤギ（はずれを意味する）がいる。プレイヤーは新車のドアを当てると新車がもらえる。プレイヤーが1つのドアを選択した後、モンティが残りのドアのうちヤギがいるドアを開けてヤギを見せる。

ここでプレイヤーは最初に選んだドアを、残っている開けられていないドアに変更しても良いと言われる。プレイヤーはドアを変更すべきだろうか？」

─⁠─「モンティ・ホール問題」『⁠Wikipedia日本語版』
2011年1月24日 12:12（UTC）

直感的には、ドアを変更してもしなくても、当たりの確率は1/2だと思えます。

この問題は非常に勘違いしやすいものであり、有名な数学者でも間違えたという逸話が残っています。こういう問題を間違えずに解くためには、理屈を考えるだけでは不十分で、乱数を使ったシミュレーションで確認するのが確実だと思われます。

この問題を忠実に再現するRubyプログラムを書くと次のようになります。

trials = (ARGV.shift || 100000).to_i # 試行回数
fail = 0 # 当たりを引いた回数
success = 0 # ハズレを引いた回数

trials.times {
  # 0、1、2の3つのドアから当たりドアをランダムに決める
  # プレイヤーは0番のドアをもらう
  hit = []
  hit[rand(3)] = true

  # モンティ・ホールの操作
  # ランダムにドア1、2を開け、ハズレが出るまで繰り返す
  while true do
    opendoor = rand(2) + 1 # ランダムに1か2を選択
    break if ! hit[opendoor]
  end
  print "モンティは#{opendoor}番のドアを開けました。"

  if hit[0] then #
    success += 1; puts "当たりです!"
  else
    fail += 1; puts "ハズレです!"
  end
}

puts "Success: #{success}"
puts "Fail: #{fail}"
puts "Ratio: #{success.to_f/trials}"

これを実行すると、選択を変更しない場合は当たりを引く確率がおよそ1/3である（選択を変更すれば当たりを引く確率が2/3に高まる）ことが明らかになります。

% ruby monty.rb
モンティは2番のドアを開けました。ハズレです!
モンティは2番のドアを開けました。ハズレです!
モンティは2番のドアを開けました。ハズレです!
モンティは1番のドアを開けました。当たりです!
モンティは1番のドアを開けました。ハズレです!
....
モンティは1番のドアを開けました。ハズレです!
モンティは1番のドアを開けました。当たりです!
Success: 3354
Fail: 6646
Ratio: 0.3354

プログラムを見ると、モンティの行為はプレイヤーの当たりハズレに何の効果も与えていないことが明白なので、モンティが何をしても当たりの確率が1/3のまま変化しないことも明らかですが、頭の中だけで考えると確率が1/2になったように勘違いしてしまうことが多いようです。

確率の計算をするとき、少しでも根拠に疑問がある場合は、乱数を用いたシミュレーションによって確認することが有効でしょう。

乱数とランダム感

「ランダムなもの」と「ランダムに感じられるもの」では微妙な違いがあります。

ランダムなのにランダムではなく感じる

完全にランダムな数値列であっても、パターンがあるように感じてしまう場合があります。

円周率

たとえば円周率の数字列は乱数列であるはずですが、小数点以下30桁以上も「0」が出現しませんし、次のような繰り返しパターンが見られるので、「⁠本当にランダムなのか？」と疑ってしまうかもしれません。

3.1415926535897932384626

ランダムな数字

次のRubyプログラムを実行して200個の数字をランダムに出力します。

200.times {
  print rand(10)
}

すると、たとえば次のような数字列が得られます。同じ数字が続いたり集中的に出現したりしている場所が意外と多いことがわかります。

5537341995738855180096821833451569453844512654843598428
0087032949714244066204568524660685406218469646838262876
9542110440797093575480604480484163727741078660949341570
92911640557823043091823428633318910

ランダムな二次元座標

また、次のようなプログラムでランダムに二次元座標を生成してプロットすると図2のようになりますが、座標の分布には偏りがあるように感じてしまいます。

1000.times {
  puts "#{rand(256)} #{rand(256)}"
}

ランダムなスライドショー

10枚の写真をランダムにスライドショー表示しようとする場合、現在の写真と同じ写真が次も表示される確率は1/10ですが、同じ写真が続けて表示されることが多いと、ランダム性が低いように感じられてしまいます。ランダムに写真が表示されるように感じるためには、同じような写真が続けて表示されないような工夫が必要です。iPodのような音楽再生機器には「シャッフル再生」という機能がありますが、本当の乱数を使ってシャッフル再生を行うと同じ曲が連続して再生されることが多々あることになり、選曲がランダムでないと感じられてしまうので、そうならない工夫がされているようです。本当にランダムな値は一様に分布することはなく、片寄りが見られる場合も多いので、人間にとってランダムに感じられるようにするためには本当の乱数を利用するのではなく、ランダムに感じられるような数字列を使うのが効果的です。

N枚のスライドをランダムに再生するとき、すべてのスライドが表示されるにはかなり時間がかかるのが普通です。何枚表示したときすべてのスライドが表示されるのかを調べるRubyプログラムは次のようになります。

n = (ARGV.shift || 100).to_i # スライドの数
trials = (ARGV.shift || 10000).to_i # 試行回数
histdiv = (ARGV.shift || 20).to_i # ヒストグラムの粒度
lim = (ARGV.shift || 2000).to_i # 表示上限
hist = [] # ヒストグラム

trials.times {
  shown = []
  nshown = 0
  count = 0
  while nshown < n do
    ind = rand(n)
    if !shown[ind] then
      shown[ind] = true
      nshown += 1
    end
    count += 1
  end
  hist[count/histdiv] = hist[count/histdiv].to_i + 1
}

100000.times { |i|
  puts "#{i*histdiv} #{hist[i].to_i}"
  break if i*histdiv >= lim
}

100枚の写真のスライドショーを行う場合、500枚ぐらいランダムに表示しても、なかなかすべての写真は表示されないことがわかります（図3⁠）⁠。

ニセ乱数

本当にランダムなものではなくランダムっぽく感じられる値を得るためには、最近出ていない値が出やすいようなニセ乱数を作るとよいでしょう。たとえば次のような工夫をします。

最近出た値は出さない
最近一度も出ていない値は出現確率を上げる

このような性質を持つニセ乱数関数rand2()を用いて、先ほどのスライドショープログラムを書き換えます。

n = (ARGV.shift || 100).to_i # スライドの数
trials = (ARGV.shift || 10000).to_i # 試行回数
histdiv = (ARGV.shift || 10).to_i # ヒストグラムの粒度
lim = (ARGV.shift || 1000).to_i # 表示上限
hist = [] # ヒストグラム

$randhist = [] # ニセ乱数の過去値リスト

def rand2(n)
  # 出やすさを決めてprob[]に格納
  prob = []
  $randhist.each { |val|
    prob[val] = 1
  }
  # 特に最近出たものは外す
  recentlen = 3
  recentlen = n-1 if recentlen >= n
  (1..recentlen).each { |i|
    break if $randhist[-i].nil?
    prob[$randhist[-i]] = 0
  }
  # 出ていないもの優先で割り当て
  (0...n).each { |i|
    if !prob[i] then
      prob[i] = 10 # 10倍出やすくする
    end
  }

  # prob[] に基づいてニセ乱数を計算
  a = []
  (0...n).each { |i|
    prob[i].times {
      a.push(i)
    }
  }
  v = a[rand(a.length)]

  # ヒストリ保存
  $randhist.shift if $randhist.length >= n
  $randhist.push(v)
  return v
end

trials.times {
  shown = []
  nshown = 0
  count = 0
  while nshown < n do
    ind = rand2(n)
    if !shown[ind] then
      shown[ind] = true
      nshown += 1
    end
    count += 1
  end
  hist[count/histdiv] = hist[count/histdiv].to_i + 1
}

100000.times { |i|
  puts "#{i*histdiv} #{hist[i].to_i}"
  break if i*histdiv >= lim
}

rand()の代わりにrand2()を使うと、図4のようなヒストグラムが得られます。160回あたりにピークがあり、300回すればほぼ確実にすべてのスライドを表示できることがわかります。

二次元ニセ乱数

ランダムに感じる二次元表示を行いたい場合も、本当の乱数を使うよりもニセ乱数を利用するほうがランダム感が出ます。

図5の「文字列探しパズル」では普通の乱数を使って文字を並べているため、同じ文字が3個以上縦や横に並んでいる場所がいくつも存在し、ランダム感があまり感じられません。

一方、前述のニセ乱数と同じような方法を使って上下と同じ文字が出にくいようにすると、図6のようになります。図5よりもランダムさが大きいように見えるでしょう。

乱数の精度

乱数で簡単にいろいろなシミュレーションができることはわかりましたが、その結果はどの程度信用できるものでしょうか？

最近のプログラミング言語には優秀な乱数生成ライブラリがあるのが普通です。たとえばRuby 1.8以降のrand()関数はメルセンヌツイスターという優秀な乱数生成手法を利用していますが、古いライブラリでは問題が出ることもあります。たとえばCの古いrand()関数でランダムな座標を生成しようとすると、条件によっては乱数として利用できません。

#include
main()
{
  for(;;){
    int x = rand();
    int y = rand();
    if(x < 0x100000 && y < 0x100000){
      printf("%d %d\n",x & 0xff,y & 0xff);
    }
  }
}

上のプログラムをSnowLeopard上で実行して生成した座標をプロットすると、図7のように平行線上に点が並んでしまいます。

ここではrand()が小さな値を返したときだけプロットするという変な条件を付けているために目に見える形でこのような不具合が発生しています。今回の記事で紹介したような簡単なシミュレーションではこういう問題が出ることはめったにありませんが、複雑なシミュレーションを行うときは、予想外の問題が発生するのを避けるため、性質の良い乱数を使うように注意するのがよいでしょう。

この問題の場合はrand()をrandom()に変更するとこのようなことは起こりません。乱数を使って複雑なシミュレーションを行う場合は注意しましょう。

連載のおわりに

連載「データ発見隊」は今回が最終回になります。1年にわたって、ちょっと変わったデータの扱いに関する話題をいろいろ紹介してきました。ふだん普通に使っているファイルやデータでも、見方を変えることによって新たな活用ができる機会はまだまだ多いと思いますので、いろんな視点でデータを活用していただければと願っています。