---
title: "Introduction to Monte Carlo simulation"
output:
  pdf_document: default
---

## Introduction to R

R is a programming language. There are often multiple ways to accomplish a given
task, so you may see different code here than in the textbook. The textbook code
is a good reference, but during the semester, don't forget your class notes and 
examples too!

Note that R is case-sensitive. That means capitalization and spelling count! The
`#` in a line of code is a commenting convention. R ignores anything after the 
`#` as a comment, but it will still show up in the line of code. 

There will be R problems in homework, plus occasionally regular textbook 
homework problems will require some R calculations (instead of using a 
calculator or tedious calculation).

## Basic built-in functions

```{r}
x <- c(3,5,6,8) # creates a vector with 4 entries `c` stands for concatenate
print(x)
sum(x) # adds up the entries in x
sum(x)/length(x) # finds the average by summing and dividing by the number of entries
mean(x) # finds the average value (just like the previous line of code)
x == 6 # == means "is"
```

## Functions useful for simulation

One of the primary reasons we will be using R and RStudio is to run simulations 
that illustrate and reinforce probability concepts. Later, it will also be 
useful for working with probability distributions. 

### The sample command

Suppose we want to generate a random number between 0 and 10 (integers only). 
There are many ways to do this, but we'll start with learning the sample 
command.

```{r}
0:5
sample(0:5, size = 1)
```


#### Question

How would you adjust the `sample` command to get 3 random values instead 
of one? What is the meaning of the `replace` option?

```{r}
## write code here
```

#### Rolling a die

To simulate rolling a standard 6-sided die:

```{r}
sample(1:6, size = 1) # roll once
sample(1:6, size = 2, replace = TRUE) # roll twice
sum(sample(1:6, 2, replace = TRUE)) # roll twice and sum the results
```

#### Tossing a coin

To simulate tossing 20 fair coins:

```{r}
sample(0:1, size = 20, replace = TRUE) # 1 corresponds to heads, 0 to tails
```

To simulate tossing 20 fair coins and counting the number of heads:

```{r}
sum(sample(0:1, size = 20, replace = TRUE))
```

### The replicate function

The replicate function lets you repeat a command a fixed number of times.

```{r}
replicate(5, print("Hello"))
```


#### Question

Write some code to repeat the experiment the following experiment 100 times:
toss a coin 20 times and count the number of heads.

```{r}
## write code here
```

### Defining a function (to make an indicator variable)

When doing Monte Carlo simulation, we want to be able to simulate an experiment
and check whether our desired event $A$ occurred. It's useful to make a function
that outputs 1 when a trial of our experiment results in event $A$ occurring,
and 0 when it does not.

The following function runs one trial of the experiment where we toss a coin 6
times and checks whether we got 3 heads.

```{r}
trial.simulation = function() { 
  # this function is like our X_k random variables
  x = sum(sample(0:1, size = 6, replace = TRUE))
  if (x == 3)
    return(1)
  else
    return(0)
} 
# when we get 3 heads in 6 tosses, this function returns 1;
# otherwise it returns 0
```

### Doing Monte Carlo simulation

Once we've written our `trial.simulation` function, to do Monte Carlo
simulation we simply need to repeat our experiment a large number of times
and average the results.

```{r}
mean(replicate(n = 1e6, trial.simulation())) 
```

#### Question
How does this compare with the exact value (found using counting)?

### Setting Seeds

Let's generate two more random numbers - just from 1 to 10, using separate commands.

```{r}
sample(1:10, 1)
sample(1:10, 1)
```

You likely got different results. Generating random numbers is great except that
we like to be able to reproduce our work. If we each use a `sample` command, we 
aren't assured of getting the same random number, unless we set the same *seed* 
for the random number generator. If you run the next code chunk (all of it), 
your randomly generated value should be the same always. 

```{r}
set.seed(1)
sample(1:10, 1, replace = TRUE)
set.seed(1)
sample(1:10, 1, replace = TRUE)
```

It may sound odd to be able to generate the same random numbers every time, but 
this is a way to make your work reproducible. You won't always need to set a 
seed, but if you find you are getting unexplainable results, and want to ask for
help, try comparing with a friend who has set the same seed.

## Getting help

What if you don't know what a command does, and there is no instruction about 
it? Help for R functions can be accessed by typing `?functionname` at the 
console. For example, if I wanted help for the sample function, I would type 
`?sample` in the console.

# Exercises

## Problem 1

Let's write code for a Monte Carlo simulation that estimates the probability of 
getting exactly 1 four in 5 rolls of a die?

#### Question
What is the exact probability (found using counting)? This will help us check 
whether our simulation is correct.

### Monte Carlo simulation

Write a trial simulation function and then use the `replicate` and `mean`
functions to get an approximation of the desired probability. Try doing
100,000 trials of the simulation. Compare with the exact probability you found.

```{r}
# write code here
```


## Problem 2

Suppose we roll a pair of dice and sum the results. What is the probability of 
getting a sum of at least 9? Write a Monte Carlo simulation to verify your 
answer.

```{r}
# write code here
```



