#-------------------------------------------------------------------------------
#                           Simulation
#-------------------------------------------------------------------------------

# Estimated parameters of the exponential distribution
x.rate <- length(data$num)

# Remember that mean = 1/x.rate
# meaning that, on average, we expect a new arrival every 1/74 of an hour.
# (1/74 =~ 0.01355)
# The random number generated using rexp() will be fractions of an hour.
lambda <- x.rate

# Number of simulated samples
n <- 74

# Set.seed()
set.seed(300)

# Let's generate random samples
# The random number generated using rexp() will be fractions of an hour.
simulated <- rexp(n,rate=lambda)

# We expect mean(simulated) =~ 0.01355* of an hour
avg <- mean(simulated)
print(paste('Average inter-arrival time (fraction of an hour)',avg)) # *
print(paste('Average inter-arrival time (minutes)',avg*60))
print(paste('Average inter-arrival time (seconds)',avg*60*60))

# Let's convert the simulated interarrival times into minutes
simulated.min <- simulated * 60

# Let's plot the histogram
breaks.points <- seq(0,max(simulated.min)+1,bins.range)
hist(simulated.min,col='chartreuse',breaks=breaks.points,border='seagreen',density=30)

boxplot(simulated.min,col='chartreuse',border='seagreen',horizontal=T,xlab='Minutes',
        main='Simulated interarrival times')

# Let's simulate the number of people per arrival
people.per.arr <- sample(x=n.of.people,size=n,replace=TRUE,prob=p.vec)

hist(people.per.arr,xlab='No. of people',breaks=break.points2,
     main='Simulated number of people per arrival',col='chartreuse',
     border='seagreen',density=30,angle=30)