```
rNx <- function(k = 1E3, x = 1) {
out <- numeric()
for (i in 1:k) {
nx <- 1
u <- runif(1)
while(u <= x) {
u <- u + runif(1)
nx <- nx + 1
}
out <- append(out, nx)
}
return(out)
}
N_1 <- rNx() # random sample of size 1000 from the distribution of N_1
mean(N_1) # estimate of 'e' by Monte Carlo
exp(1) # value of 'e'
```