```
<- function(seed, n, min = 0, max = 1) {
runif_zx81 <- 2^16 + 1
m if(seed > m) return('Error: seed cannot exceed 65537.')
if(seed <= 0) return('Error: must be positive.')
<- 75
a <- vector(mode = 'numeric', length = n+1)
x 1] <- seed
x[for(i in 1:n) {
+ 1] <- (a * x[i]) %% m
x[i
}+ (x[-1] / m) * (max - min)
min }
```

# Exercise - Monte Carlo On a Desert Island

When it comes to basic numerical methods, “rolling your own” is generally a bad idea. If you need to generate (pseudo)random numbers in everyday life, you should rely on the excellent work of generations of numerical analysts and computer scientists and just use the base R functions `runif()`

, `rnorm()`

etc. But what if you were stranded on a desert island and needed to re-create civilization from scratch? Surely you would need a way to simulate uniform and normal random variables! In this exercise you will generate pseudorandom numbers as though you were Robinson Crusoe. Again: unless civilization collapses, you *you should not do this in real life*. Nevertheless, this exercise will give you a better appreciation for how much work goes on behind the scenes when you run a simulation study!

- The ZX81 was a personal computer sold in the UK in the early 1980s. To simulate standard uniform draws, the ZX81 relied on the Lehmer random number generator discussed in our lecture on Monte Carlo simulation. In particular, it set \(m = 2^{16} + 1\) and \(a = 75\).
- Write a function called
`runif_zx81()`

that generates uniform draws in the same way as the Z81. Your function should take four arguments. The first argument,`seed`

, should be the seed of the Lehmer random number generator. The remaining arguments–`n`

,`min`

and`max`

–should be identical to the corresponding arguments of the base R function`runif()`

. - Choose a seed and make 1000 random draws from
`runif_zx81()`

. Do these appear to be iid Uniform(0,1) draws? Make a histogram, QQ plot, and time series plot to check.

- Write a function called
- Let \(U_1\) and \(U_2\) be a pair of independent standard uniform random variables. Now define two new random variables: \(R \equiv \sqrt{-2 \log(U_1)}\) and \(\Theta \equiv 2\pi U_2\). Then it can be shown that \(Z_1 \equiv R \cos(\Theta)\) and \(Z_2 \equiv R \sin(\Theta)\) are a pair of independent standard normal random variables. This result is usually called the Box-Muller Transform.
- Write a function called
`rnorm_zx81()`

that uses`runif_zx81()`

along with the the Box-Muller transform to simulate draws from a normal random variable. The first argument,`seed`

, should be the seed that you pass to`runif_zx81()`

. The remaining arguments–`n`

,`mean`

, and`sd`

–should match those of the base R function`rnorm()`

. - Choose a seed and make 1000 random draws from
`rnorm_zx81()`

. Do these appear to be iid standard normal draws? Make a histogram, QQ plot, and time series plot to check.

- Write a function called

## Solutions

### 1(a)

### 1(b)

```
library(tidyverse)
<- runif_zx81(1983, 1000)
u tibble(u) |>
ggplot(aes(u)) +
geom_histogram(bins = 20)
```

```
tibble(u) |>
ggplot(aes(1:length(u), u)) +
geom_line()
```