# Ambulances


Find this notebook on the web at
<a class="quarto-xref" href="https://resampling-stats.github.io/latest-r/resampling_method.html#nte-ambulances">Note <span>2.1</span></a>.

**Note: Comments in code**

The first thing to say about the code you will see below is there are
some lines that do not do anything; these are the lines beginning with a
`#` character (read `#` as “hash”). Lines beginning with `#` are called
*comments*. When R sees a `#` at the start of a line, it ignores
everything else on that line, and skips to the next. Here’s an example
of a comment:

```{r}
# R will completely ignore this text.
```

Because R ignores lines beginning with `#`, the text after the `#` is
just for us, the humans reading the code. The person writing the code
will often use comments to explain what the code is doing.

**End of Note: Comments in code**

Our next task is to use R to simulate a single day of ambulances. We
will again represent each ambulance by a random number from 0 through 9.
16 of these numbers represents a simulation of all 16 ambulances
available to the contractor. We call a simulation of all ambulances for
a specific day — one *trial*.

Recall that we want 16 10-sided dice — one per ambulance. Our dice
should be 10-sided, because each ambulance has a 1-in-10 chance of being
out of order.

The program to simulate one trial of the ambulances problem therefore
begins with these commands:

```{r}
# Ask R to generate 16 numbers from 0 through 9.

# These are the numbers we will ask R to select from.
numbers <- c(0, 1, 2, 3, 4, 5, 6, 7, 8, 9)

# Get 16 values from the *numbers* sequence.
# Store the 16 numbers with the name "a"
# We will explain the replace=TRUE later.
a <- sample(numbers, 16, replace=TRUE)

# The result is a sequence of 16 numbers.
a
```

The commands above ask the computer to store the results of the random
drawing in a location in the computer’s memory to which we give a name
such as `a` or `ambulances` or `aardvark` — the name is up to us.

Next, we need to count the number of defective ambulances:

```{r}
# Count the number of nines in the random numbers.
# The "a == 9" part identifies all the numbers equal to 9.
# The "sum" part counts how many numbers "a == 9" found.
b <- sum(a == 9)
# Show the result
b
```

You see that we have code like this:

```{r}
sum(a == 9)
```

to count the number of values equal to 9 in the sequence `a`. For now
just read sum(a == 9 as “count the number of 9s in `a`”. We will explain
how this code works fairly soon
(<a class="quarto-xref" href="https://resampling-stats.github.io/latest-r/resampling_with_code.html#sec-counting-results"><span>Section 5.12</span></a>).

The `sum` command is a *counting* operation. It asks the computer to
*count* the number of `9`s among the 16 numbers that are in location `a`
following the random draw carried out by the `sample` operation. The
result of the `sum` operation will be somewhere between 0 and 16, the
number of simulated ambulances that were out-of-order on a given
simulated day. The result is then placed in another location in the
computer’s memory that we label `b`.

Above you see that we have worked out how to tell the computer to do a
single trial — one simulated day.

### 2.3.1 Repeating trials

We could run the code above for one trial over and over, and write down
the result on a piece of paper. If we did this 100 times we would have
100 counts of the number of simulated ambulances that had broken down
for each simulated day. To answer our question, we will then count the
number of times the count was three or more, and divide by 100, to get
an estimate of the proportion of days with three or more out-of-order
ambulances.

One of the most useful things about the computer is that it is very good
at repeating tasks many times. Our next job is to ask the computer to
repeat the single trial many times — say 1000 times — and count up the
results for us.

Of course R is very good at repeating things, but the instructions to
tell R to repeat things will take a little while to get used to. Soon,
we will spend some time going over it in more detail. For now though, we
show you how what it looks like, and ask you to take our word for it.

The standard way to repeat steps in R is a `for` loop. For example, let
us say we wanted to display “Hello” five times. Here is how we would do
that with a `for` loop:

```{r}
# Read the next line as "repeat the following steps five times".
for (i in 1:5) {
    # The stuff between the curly brackets is the code we
    # repeat five times.
    # Print "Hello" to the screen.
    message("Hello")
}
```

You can probably see where we are going here. We are going to put the
code for one trial inside a `for` loop, to repeat that trial code many
times.

Our next job is to *store* the results of each trial. If we are going to
run 1000 trials, we need to store 1000 results.

To do this, we start with a sequence of 1000 zeros, that we will fill in
later, like this:

```{r}
# Ask R to make a sequence of 1000 zeros that we will use
# to store the results of our 1000 trials.
# Call this sequence "z"
z <- numeric(1000)
```

For now, `z` contains 1000 zeros, but we will soon use a `for` loop to
execute 1000 trials. For each trial we will calculate our result (the
number of broken-down ambulances), and we will store the result in the
`z` store. We end up with 1000 trial results stored in `z`.

With these parts, we are now ready to solve the ambulance problem, using
R.

### 2.3.2 The solution

This is our big moment! We will combine the elements shown above to
perform our ambulance simulation over, say, 1000 days. Just a quick
reminder: we do not expect you to understand all the detail of the code
below; we will cover that later. For now, see if you can follow along
with the gist of it.

To solve resampling problems, we typically proceed as we have done
above. We figure out the structure of a single trial and then place that
trial in a `for` loop that executes it multiple times (once for each
day, in our case).

Now, let us apply this procedure to our ambulance problem. We simulate
1000 days. You will see that we have just taken the parts above, and put
them together. The only new part here, is the step at the end, where we
store the result of the trial. Bear with us for that; we will come to it
soon.

```{r}
# Ask R to make a sequence of 1000 zeros that we will use
# to store the results of our 1000 trials.
# Call this sequence "z"
z <- numeric(1000)

# These are the numbers we will ask R to select from.
numbers <- 0:9

# Read the next line as "repeat the following steps 1000 times".
for (i in 1:1000) {
    # The stuff between the curly brackets is the code we
    # repeat 1000 times.

    # Get 16 values from the *numbers* sequence.
    # Store the 16 numbers with the name "a"
    a <- sample(numbers, 16, replace=TRUE)

    # Count the number of nines in the random numbers.
    # The "a == 9" part identifies all the numbers equal to 9.
    # The "sum" part counts how many numbers "a == 9" found.
    b <- sum(a == 9)

    # Store the result from this trial in the sequence "z"
    z[i] <- b

    # Now go back and repeat the trial, until done.
}
```

The `z[i] &lt;- b` statement that follows the `sum` *counting* operation
simply keeps track of the results of each trial, placing the number of
defective ambulances for each trial inside the sequence called `z`. The
sequence has 1000 positions: one for each trial.

When we have run the code above, we have stored 1000 trial results in
the sequence `z`. These are 1000 counts of out-of-order ambulances, one
for each of our simulated days. Our last task is to calculate the
proportion of these days for which we had more than three broken-down
ambulances.

Since our aim is to count the number of days in which three or more
defective ambulances occur, we use another *counting* `sum` command at
the end of the 1000 trials. This command *counts* how many times three
or more defects occurred in the 1000 days recorded in our `z` sequence,
and we place the result in another location, `k`. The value of `k` gives
us the total number of days where 3 or more defective ambulances are
seen to occur. Then we divide the number in `k` by 1000, the number of
trials. Thus we obtain an estimate of the chance, expressed as a
probability between 0 and 1, that 3 or more ambulances will be defective
on a given day. And we store that result in a location that we call
`kk`, which R subsequently prints to the screen.

```{r}
# How many trials resulted in 3 or more ambulances being out of order?
# (The test counts the number of days where the out-of-order count was greater
# than or equal to 3, written as ">= 3").
k <- sum(z >= 3)

# Convert to a proportion.
kk <- k / 1000

# Show the result.
message(kk)
```

This is the estimate we wanted; the proportion of days where three or
more ambulances were out of action.

We have crept up on the solution, so it might not be clear to you how
few steps you needed to do this task. Here is the whole solution to the
problem, without the comments:

```{r}
z <- numeric(1000)
numbers <- 0:9

for (i in 1:1000) {
    a <- sample(numbers, 16, replace=TRUE)
    b <- sum(a == 9)
    z[i] <- b
}

k <- sum(z >= 3)
kk <- k / 1000
message(kk)
```
