\documentclass{article}

\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{graphics}
\usepackage{listings}
\usepackage{tcolorbox}
\usepackage{verbatim}

\DeclareMathOperator*{\argmin}{arg\,min}

\marginparwidth 0pt
\oddsidemargin  0pt
\evensidemargin  0pt
\marginparsep 0pt
\topmargin   0pt
\textwidth   6.5in
\textheight  8.5in
\setlength{\parindent}{0em}

\newenvironment{packed_itemize}{
\begin{itemize}
  \setlength{\itemsep}{0pt}
  \setlength{\parskip}{0pt}
  \setlength{\parsep}{0pt}
}{\end{itemize}}

\newenvironment{packed_enumerate}{
\begin{enumerate}
  \setlength{\itemsep}{1pt}
  \setlength{\parskip}{0pt}
  \setlength{\parsep}{0pt}
}{\end{enumerate}}

\definecolor{codeback}{RGB}{235,255,255}

\lstset{ %
language=Python,                % choose the language of the code
basicstyle=\footnotesize,       % the size of the fonts that are used for the code
numbers=left,                   % where to put the line-numbers
numberstyle=\footnotesize,      % the size of the fonts that are used for the line-numbers
stepnumber=1,                   % the step between two line-numbers. If it is 1 each line will be numbered
numbersep=10pt,                  % how far the line-numbers are from the code
backgroundcolor=\color{codeback},  % choose the background color. You must add \usepackage{color}
showspaces=false,               % show spaces adding particular underscores
showstringspaces=false,         % underline spaces within strings
showtabs=false,                 % show tabs within strings adding particular underscores
frame=single,           % adds a frame around the code
tabsize=2,          % sets default tabsize to 2 spaces
captionpos=b,           % sets the caption-position to bottom
breaklines=true,        % sets automatic line breaking
breakatwhitespace=false,    % sets if automatic breaks should only happen at whitespace
escapeinside={(*@}{@*)}          % if you want to add a comment within your code
}

\begin{document}

\begin{center}
  \Huge{Monte Carlo Methods}\\
  \Large{MATH1900: Machine Learning}
\end{center}
Location: http://people.sc.fsu.edu/$\sim$jburkardt/classes/ml\_2019/monte\_carlo/monte\_carlo.pdf
\hrule

\begin{center}
  \includegraphics[scale=0.60]{random_selection.png}
  \vskip 0.1in
  {\it{Sampling a big data set.}}
\end{center}
\vskip 0.1in

A Monte Carlo method is a computational process that uses random sampling.  This lab
introduces a few examples, in which we estimate averages, areas, and integrals.  

\section{Example: Estimating statistics of a very large set}

A surveys claims that 68\% of the population are fans of the Peppa Pig show, and
that each such viewer watches for an average of 23 minutes a day.  These claims 
were not made by asking every person in the country, but rather, by selecting a 
random sample of the population, small enough to treat individually, but large 
enough to be representative.  In machine learning, it may be the case that a 
dataset contains millions of entries, and that each query is rather expensive.  
If we only need a rough estimate for our result (say, just two digits of
accuracy), then random sampling may give us an answer that is good enough.
\vskip 0.1in
For simplicity, we will assume that our data is simply a vector $x$ of $n$ real 
values, and that we want the minimum, average, and maximum.  Instead of 
checking all $n$ values, we will choose a much smaller number of values $k$.  
To carry out our task, we need to generate $k$ distinct random indices $i$ 
such that $0 \leq i < n$.  We do this with {\tt{np.random.choice()}}:
\begin{lstlisting}
  import numpy as np
  n = 1000
  mu = 23
  sigma = 4
  x = mu + sigma * np.random.randn ( n )
  x_min = np.min ( x )
  x_mean = np.mean ( x )
  x_max = np.max ( x )
  
  k = 5
  i = np.random.choice ( n, k, replace = False )
  xs = x[i]
  xs_min = np.min ( xs )
  xs_mean = np.mean ( xs )
  xs_max = np.max ( xs )
\end{lstlisting}

\section{Exercise: Estimating height and weight statistics}

The file {\it{hw\_25000.txt}} contains 25,000 records from the physical exams
of a population.  Each record consists of three values: an index (from 1 to 25,000), 
a height (inches), and a weight (pounds).  Suppose we needed to estimate the
minimum, mean, maximum and standard deviation of the height and weight of the
entire population, but cannot afford to view all the records.  Estimate
these 4 statistical quantities for height and for weight, using random
samples of size 10, 100, and 1,000.  Compare with the exact values.
\begin{verbatim}
                   <-----Height----->                          <-----Weight----->
     N     Min       Mean       Max        Std         Min       Mean       Max        Std
        ---------  ---------  ---------  ---------  ---------  ---------  ---------  ---------
    10  .........  .........  .........  .........  .........  .........  .........  .........
   100  .........  .........  .........  .........  .........  .........  .........  .........
 1,000  .........  .........  .........  .........  .........  .........  .........  .........
25,000  .........  .........  .........  .........  .........  .........  .........  .........
\end{verbatim}

\section{Example: Estimate the integral of the ``humps'' function}

The mean value theorem for integrals says
\begin{displaymath}
  \int_a^b f(x) \, dx = \bar{f} * ( b - a )
\end{displaymath}
where $\bar{f}$ is the average value of $f(x)$ over $[a,b]$.  This gives us a way
to estimate an integral by sampling.  Simply choose a random sample of points $x$ in
$[a,b]$, average the values $f(x)$, and multiply by $(b-a)$.  As the number of points
used increases, we improve our estimate of $\bar{f}$ and hence of the integral.
\begin{lstlisting}
  def humps ( x )
    y = 1.0 / ( ( x - 0.3 )**2 + 0.01 ) + 1.0 / ( ( x - 0.9 )**2 + 0.04 ) - 6.0
    return y
    
  import numpy as np
  
  a = 0.0
  b = 2.0
  exact = 29.32621380439115
  
  for logn in range ( 2, 6 ):
    n = 10**logn
    x = a + ( b - a ) * np.random.rand ( n )
    fx = humps ( x )
    estimate = ( b - a ) * np.mean ( fx )
    err = exact - estimate
    print ( 'n =', n, 'estimate =', estimate, 'error =', err )
\end{lstlisting}

\section{Exercise: Estimate the integral of $\frac{1}{\sqrt{x}}$}

What is the exact value of $\int_0^4 \frac{1}{\sqrt{x}} dx$?
\vskip 0.1in
Use a Monte Carlo method to estimate this quantity.  Choose a value of $n$ large
enough that your error is less than 0.001.
\begin{verbatim}
Exact integral value = ..........
Estimated value      = ..........  using N = .......... samples
Error                = ..........
\end{verbatim}

\section{Example: Estimating the area of an ellipse}

If an ellipse has the formula:
\begin{displaymath}
\frac{x^2}{a^2} + \frac{y^2}{b^2} = 1- 
\end{displaymath}
then it lies inside the rectangle $-a \leq x \leq a, -b \leq y \leq b$.  Hence,
the area inside the ellipse must be no greater than $4*a*b$.  How could we estimate
this area?  If we randomly sample the surrounding rectangle, then it is reasonable,
(and correct) to estimate the ellipse area by
\begin{displaymath}
\text{ellipse area} \approx \frac{\text{points inside ellipse}}{\text{points inside rectangle}} * \text{area of rectangle}
\end{displaymath}
So we can generate $n$ points inside the rectangle, and then we can determine if each
 point is inside the ellipse by checking if $\frac{x^2}{a^2} + \frac{y^2}{b^2} \leq 1$.
\begin{lstlisting}
import numpy as np
a = 3
b = 2

n = 100
x = - a + 2 * a * np.random.rand ( n )
y = - b + 2 * b * np.random.rand ( n )
i = np.count_nonzero ( x**2/a**2 + y**2/b**2 <= 1 )
area_estimate = i / n * 4 * a * b
\end{lstlisting} 
 
It turns out that this area is exactly $pi*a*b$.  Can you modify this procedure 
so that it increases $n$ until the error is less than 0.01?

\section{Exercise: Estimating the area of a circle minus a circle}

Let $C_1$ be the unit circle of radius 1 and center (0.0,0.0), and $C_2$ be the
circle of radius 0.1 and center (0.5,0.0).  We want to estimate the area
of $C_3$, the shape formed by $C_1$ ``minus'' $C_2$, in other words, $C_1$ with
a hole in it.
\vskip 0.1in
\begin{center}
  \includegraphics[scale=0.30]{circle_hole.png}
  \vskip 0.1in
  {\it{The $C_3$ region.}}
\end{center}
\vskip 0.1in
What are the areas of $C_1$, $C_2$, and $C_3$?
\vskip 0.1in
Use a Monte Carlo method to estimate the area of $C_3$.  Try to use a value of
$n$ large enough that your error is less than 0.01.  Note that, in your
call to {\tt{np.count\_nonzero ( )}} you now need {\bf{two}} conditions to be true
in order that a point is inside $C_1$ and not inside $C_2$.
You specify these, using an ampersand, as in:
\begin{lstlisting}
i = np.count_nonzero {\tt{( condition 1 & condition 2 )}}.
\end{lstlisting}
Think for a moment about what condition 1 and condition 2 should be!
\vskip 0.1in
\begin{verbatim}
    exact C_1 area .........
    exact C_2 area .........
    exact C_3 area .........
estimated C_3 area .........  using N = .........
             error .........
\end{verbatim}

\end{document}
