ex1 <- function(total = 24, marked = 4, boxes = 4, iterations = 10000) {
  perbox <- total %/% boxes
  successes <- 0
  
  for (i in 1:iterations) {
    apples <- sample(append(rep("marked", each=marked), rep("normal", each=total - marked)))
    
    current <- 0
    for (b in 0:boxes-1) {
      current <- current + (1 == sum(apples[1:perbox + b*perbox] == "marked"))
    }
    
    successes <- successes + (current == 4)
  }
  
  return (round(successes / iterations, digits=2))
}

ex1()

ex2 <- function(n, lambda = 1/5, N = 10000, mu = 5, sigma = 5 / sqrt(n)) {
  x <- function() rexp(n, rate=lambda)
  xdash <- replicate(N, sum(x()) / n)

  fchapeau <- function(t) {
    # I е излишна на практика, но гледам да следвам условието възможно най-стриктно
    I <- function(a) as.integer(a == TRUE)
    return (1/N * sum(I(xdash <= t)))
  }

  t <- seq(2, 8, by=0.1)
  Fpoints <- sapply(t, fchapeau)
  Npoints <- sapply(t, function(t) pnorm(t, mu, sigma))
  
  plot (t, Fpoints, type="l")
  lines(t, Npoints, lty=2, col="darkred")
}

ex2(30)
ex2(120)
ex2(200)