(fill in your two names here)
Facilitator: (fill in name)
Spokesperson: (fill in name)
Process Analyst: (fill in name)
Quality Control: (fill in name)
If there are only three people in your team, have one person serve as both spokesperson and process analyst for the rest of this activity.
At the end of this Lesson, you will be asked to record how long each Model required for your team. The Facilitator should keep track of time for your team.
The Monte Carlo (MC) method is a way to approximate solutions to a variety of mathematical problems by conducting statistical sampling experiments on a computer. In other words, this method relies on random sampling to obtain numerical results. This technique was named after a city in Monaco where the roulette wheel was used essentially as a random number generator. You can find MC methods used in everything from economics to modeling molecules to regulating the flow of traffic. Although the way MC algorithms are applied varies widely from field to field, all "Monte Carlo" experiments use random numbers to examine some problem.
Here is an example of a Monte Carlo method that generates an approximation for $\pi$. Imagine you are throwing darts randomly at a square. If you were to inscribe a circle within the square, it is possible to calculate the value of $\pi$ based the proportion of darts that land inside the circle in comparison to the total number of darts thrown.
1. The square on the left in Model 1 has been divided into four equal sized squares.
1a. If you throw $n$ total darts randomly at this square, what is the probability that $m$ darts will land in the upper-right quadrant (where $m$ < $n$)? Explain.
Write your answer in terms of an expression with $m$ and $n$.
1b. Show how you can arrive at this relationship in terms of the ratio of the area of the small square and the large square.
1c. Generalize this relationship in terms assuming the area of the small square is $A$.
2. Now consider the square with the circle inscribed (shown on the right of Model1):
2a. What is the ratio the area of the circle divided by the area of the square? (LaTeX and code are both good here)
2b. What is the probability $m$ darts lie within the circle, if $n$ darts are thrown in total? Write your answer in terms of an expression with $m$ and $n$.
2c. Generalize this relationship (which will be used in the next question to show how the value of $\pi$ could be calculated).
2d. Explain how an approximation for the value of $\pi$ can be found in terms of darts.
3. Given that the center of the circle is the origin (i.e. (0, 0)), add a point $(x,y)$ inside the circle. The distance ($d$) from the origin of the circle to that point is the hypotenuse of the right triangle with sides of $x$, $y$, and $d$.
3a. On a piece of paper, sketch this right triangle over the diagram in Model 1.
3b. Provide an equation to calculate the distance, $d$, from the circle origin to any point $(x,y)$.
3c. Write a Boolean condition that returns True if a point $(x, y)$ is inside the circle or on its boundary.
4. Given the scenario of throwing random darts at a square, how many random numbers are necessary to represent each dart and what would each random number represent?
5. Consider what control structure(s) you would use (sequential, loop, condition) to implement this simulation and write high-level pseudocode for this Monte Carlo simulation to approximate $\pi$. Your function should take one parameter - the number of "darts". Show your answers to an instructor before proceeding, but do not spend more than five minutes on this step.
In the previous lesson, we generated a non-uniform distribution of random numbers by transformation. There are cases, however, when the desired distribution may not be known in analytic form. For example, suppose we want the Gaussian distribution:
$$P(y)dy = e^{-y^2} dy$$For this function the integral cannot easily be evaluated analytically. Such problems can be solved using an alternative algorithm known as the Rejection Method. This approach uses a similar Monte Carlo approach as the algorithm used above to estimate $\pi$.
Suppose you want to generate numbers that follow a certain distribution $P(y)$ with $y$ in the interval $[a,b)$. The easiest way to achieve this (although perhaps not the most efficient) is to throw away all events that are "outside" of $P(y)$ in a two dimensional plane (as we did for $\pi$). For example, if your distribution is $P(y) = e^{-y}$ in the range [0,10), as shown below:
For your $y$ value, you would use a uniformly distributed random number between $[0,10)$. As your $x$ value (i.e. second coordinate for your plane and each "dart" that you throw), you would use a second independent random number uniformly distributed between 0 and 1. (Note: this assumes a normalized $P(y)$.) If all random numbers were accepted, this would correspond to a plane of numbers uniformly distributed in the two dimensions. However, if we reject all $y$ numbers for which $x > P(y)$ (in our case $e^{-y}$), the remaining accepted $y$ will be distributed according to $e^{-y}$.
If we go back to our random dart logic, and we throw a bunch of random darts on the scales discussed above (represented by the grey points below) our system might look something like:
Now we are trying to program a condition that will determine which of the "dart throws" above should be "rejected" as not being under the curve, and change the colors of the dots that fall "inside" (and we want to keep (in blue)) vs. "outside" (and we want to reject (in orange)) of the distribution, we would get something like:
Then if we plot a histogram of the values that are not rejected, we get an estimate of our distribution.
6. Given the equation $P(y) = e^{-y}$, which of the following $x,y$ points (dart coordinates) should be rejected?
| dart x | dart y | P(y) | Accept or Reject? |
|---|---|---|---|
| 0.46 | 5.4 | ||
| 0.26 | 1.2 | ||
| 0.81 | 4.1 |
7. Write a Boolean condition that defines if a point is accepted as part of the probability distribution.
8. Write higher-level pseudocode to generate and return 100 random points that fall under a distribution P(y) using the rejection method.
A. Define and test a function called monte_pi that takes one integer parameter (the number of darts thrown, called $n$ above, but for extra laughs call this parameter thon, get it?!?!?!?) and returns the estimated value of $\pi$ using the Monte Carlo method (from Model 1).
In [ ]:
9. What happens as you increase the number of "darts"? Explain.
B. Define and test an exp_method2 that takes as a parameter the number of random numbers to generate. Use the rejection method (from Model 2) to obtain exponentially distributed random numbers between [0, 10).
In [ ]:
C. Copy and modify your exp_method2 code from above to make exp_method2_hist that will allow you to visualize the results of your function as a histogram to confirm the generated random numbers approximate an exponential distribution. Test this function with several test cases increasing the number of "darts".
In [ ]:
10. What happens as you increase the number of "darts"? Explain.
In [ ]:
How much time did it require for your team to complete each Model?
Model 1:
Model 2:
Programming: