Untitled

final_1 <- function(N, M, param, k, lambda, B, myfun="normal", mystat="median"){
  coverage_sim <- 0
  length_sim <- 0
  results <- c()
  for(n in N){
  for(m in M){
    
    # Generate m random Samples
    for(i in 1:m){
      #generate random sample 
      x = rweibull(n,shape=k,scale=lambda)
      #print(summary(x))
      # determine the actual population parameter - median
      boot.samples = bootsampling(x, boot.replicates = B)
      boot.out <- switch(
        mystat,
        "median"=boot.stats(boot.samples, median),
        "variance"=boot.stats(boot.samples, var)
        
      )
      sample_stat <- switch(
        mystat,
        "median" = median(x),
        "variance" = var(x)
      ) 
      confint_sim = switch(
        myfun,
        "normal" = normal_ci(sample_stat, boot.out),
        "basic" = basic_ci(sample_stat, boot.out),
        "percentile" = percentile_ci(sample_stat, boot.out),
        "studentized" = studentized_ci(boot.samples, boot.out, sample_stat)
      )
      #cat("confidence interval: [", confint_sim[1], ",", confint_sim[2], "]", "\n")
      length <- confint_sim[2] - confint_sim[1]
      length_sim <- length_sim + length
      if(param <= confint_sim[2] && param >= confint_sim[1]){
        coverage_sim <- coverage_sim +1
      }
      
    }
    coverage <- coverage_sim/m
    length <- length_sim/m
    #cat("n =", n, ", M =",m, ":", "\n")
    #cat("Avg. Length:", length, "Coverage %:", coverage, "\n")
    coverage_sim <- 0
    length_sim <- 0
    result <- c(coverage, length)
    results <- c(results, result)
  }
  }
  results
  }