reps <- 1000

bias_med <- NULL
Var_med <- NULL
MSE_med <- NULL

bias_mean <- NULL
Var_mean <- NULL
MSE_mean <- NULL

populationMean <- 1/2

for (n in c(10,20,40,80))
{
  Z <- array(runif(n*reps), c(n,reps))
  Med <- apply(Z, 2, median)
  Mean <- apply(Z, 2, mean)

  bias_med <- c(bias_med, mean(Med) - populationMean)
  Var_med <- c(Var_med, var(Med))
  MSE_med <- c(MSE_med, mean((Med-populationMean)^2))

  bias_mean <- c(bias_mean, mean(Mean) - populationMean)
  Var_mean <- c(Var_mean, var(Mean))
  MSE_mean <- c(MSE_mean, mean((Mean-populationMean)^2))
}




reps = 100
nboot = 1000

theta.POP = 2

SS = c(5,10,20)
Rel.Bias = NULL
for (k in 1:3) 
{
  n = SS[k]
  rel.bias = NULL
  for (r in 1:reps) 
  {
    ## Generate n observations from Unif(theta, 10)
    X = (10-theta.POP)*runif(n)+theta.POP
    
    ##  Get the "theta" for the RB formula
    theta.X  = min(X)
    
    ##  Generate the artificial data from X
    ii = ceiling(n*runif(nboot*n))
    Xboot = X[ii]
    Xboot = array(Xboot, c(nboot, n))
    
    ##  Compute nboot "thetaHat"'s for the RB formula
    ##  from this artificial data
    theta.BOOT = apply(Xboot, 1, min)
    
    ##  Compute the Relative Bias formula
    rel.bias[r] = (mean(theta.BOOT)-theta.X)/theta.X
  }
  Rel.Bias[k] = mean(rel.bias)
}






reps = 100
nboot = 1000

theta.POP = 2

SS = c(5,10,20)
Rel.Bias = NULL
for (k in 1:3) 
{
  n = SS[k]
  rel.bias = NULL
  for (r in 1:reps) 
  {
    ## Generate n observations from Unif(theta, 10)
    X = (10-theta.POP)*runif(n)+theta.POP
    
    ##  Get the estimate of "theta" 
    theta.X  = min(X)
    
    ##  Generate the artificial Uniform data using this estimate
    Xboot = (10-theta.X)*runif(nboot*n)+theta.X
    Xboot = array(Xboot, c(nboot, n))
    
    ##  Compute nboot "thetaHat"'s for the RB formula
    ##  from this artificial data
    theta.BOOT = apply(Xboot, 1, min)
    
    ##  Compute the Relative Bias formula
    rel.bias[r] = (mean(theta.BOOT)-theta.X)/theta.X
  }
  Rel.Bias[k] = mean(rel.bias)
}