This post originated from an RSS feed registered with Ruby Buzz by Eigen Class. Original Post: Warm fuzzy things for random simulations Feed Title: Eigenclass Feed URL: http://feeds.feedburner.com/eigenclass Feed Description: Ruby stuff --- trying to stay away from triviality. Latest Ruby Buzz Posts Latest Ruby Buzz Posts by Eigen Class Latest Posts From Eigenclass

Let's talk about random experiments. The simplest one is tossing a coin, with outcomes "heads" and "tails". It's so elementary that we fully understand it intuitively once we're told that the coin is not loaded. What's "everything", though? There are three basic problems of interest:

1. How to simulate the experiment with a computer and obtain results undistinguishable from the real thing?
2. If we're given an outcome, can we compute its probability?
3. Can we enumerate all possible outcomes and their probabilities?

In our first experiment, the answers are "rand 2" (if take e.g. "heads" = 0, "tails" = 1), "1/2" and "heads, tails, both with probability 1/2". Easy indeed. Now on to something harder: rolling a fair die. The answers are now "1 + rand(6)", "1/6", and "1..6, all with probability 1/6". Still trivial. Let's aim for the stars! What about rolling three dice? The outcomes are clearly 3..18, but what are the probabilities?

Intuition doesn't help much here --- we "can see" the possible outcomes, but most people cannot give the precise probabilities intuitively (I certainly can't; I can at most draw an approximate graph, but that's because I am used to convolutions). The first problem was how to simulate the experiment; this is by far the easiest one. This is what our dice rolling simulation looks like in Ruby:

  def roll_3d6
    d1 = rand(6)
    d2 = rand(6) 
    d3 = rand(6)
    1 + d1 + 1 + d2 + 1 + d3
  end

We can now simulate the experiment easily, but this doesn't help us with the harder problems: giving the probability of an outcome and doing so for all possible outcomes. Clearly, the above snippet has got all the information we need about this random simulation; it's indeed a precise definition. Wouldn't it be nice to have the computer solve the three problems of interest (simulation, distribution and expectation) for us if we give it the specification of the simulation? It turns out this can be done with a little code massaging:

  def roll_3d6(m)
    m.rand(6).in do |d1|
      m.rand(6).in do |d2|
        m.rand(6).in do |d3|
          m.ret(1+d1 + 1+d2 + 1+d3)
        end
      end
    end
  end

Given that code and suitable Simulation, Distribution and Expectation classes, we can solve all three problems with the same code.

Simulating the experiment:

 p roll_3d6(Simulation.init).run    
 #  >> 9

All the possible outcomes and their probabilities:

 p roll_3d6(Distribution.init).run  
 # >>  
 # [[3, 0.00462962962962963], [4, 0.0138888888888889], [5, 0.0277777777777778], 
 # [6, 0.0462962962962963], [7, 0.0694444444444444], [8, 0.0972222222222222], 
 # [9, 0.115740740740741], [10, 0.125], [11, 0.125], [12, 0.115740740740741],
 # [13, 0.0972222222222222], [14, 0.0694444444444444], [15, 0.0462962962962963],
 # [16, 0.0277777777777778], [17, 0.0138888888888889], [18, 0.00462962962962963]]

Getting the probability given an outcome:

 (1..21).each{|n| p roll_3d6(Expectation.init).run(n) }
 # >>
 # 0.0
 # 0.0
 # 0.00462962962962963
 # 0.0138888888888889
 # ...

What's really nice is that Simulation, Expectation and Distribution can be reused for other random experiments. Say you get an extra roll if any of the first 3 dice gave a 6:

  def roll_3d6_b(m)
    m.rand(6).in do |d1|
      m.rand(6).in do |d2|
        m.rand(6).in do |d3|
          if [d1, d2, d3].include?(5)
            m.rand(6).in{|d4| m.ret(1+d1 + 1+d2 + 1+d3 + 1+d4) }
          else
            m.ret(1+d1 + 1+d2 + 1+d3)
          end
        end
      end
    end
  end

Here's the distribution we get with roll_3d6_b(Distribution.init).run: