### Problem 1
### FOR Uniform(0,1) data
### Estimate the population mean using the sample median and sample mean
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))
}



### Problem 2
### FOR Uniform(0,1) data
### Estimate the population median using the sample median and sample mean
reps <- 1000

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

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

populationMedian <- 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) - populationMedian)
  Var_med <- c(Var_med, var(Med))
  MSE_med <- c(MSE_med, mean((Med-populationMedian)^2))

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


## From the lecture notes, here are some hints 
##  for modification of the code for homework #3

## Number of clusters.
nc <- 100
## Number of observations per cluster.
cs <- 5
## Shrinkage factors.
F <- c(0.5, 0.6, 0.7, 0.8, 0.9, 1)    # Modification Needed here
                                      # We want F to vary from 0 to 1 by 0.05
				     # consider using the seq() command

## Save all the estimates.
Estimate <- array(0, c(1000,6))       # Modification
					# Want to save an estimate for
				 	# for each cluster, for each lambda,
				      # and each repetition
				      # Consider adding a dimension to the array
				     # that corresponds to each cluster
## True cluster means.
CM <- rnorm(nc)

## Simulation replications.
for (r in 1:1000)
{
  ## Estimated cluster means.
  Z <- array(0, nc)
  for (k in 1:nc)
  {
    z <- CM[k] + rnorm(cs)
    Z[k] <- mean(z)
  }

  ## The ’grand mean’
  GM <- mean(Z)
  for (k in 1:length(F))
  {
    ## Modify the cluster 1 mean estimate by shrinking toward   ## change this
    ## the grand mean.
    S <- GM + F[k]*(Z[1]-GM)         ## Modification Needed here
    Estimate[r,k] <- S			 ## Modification Needed here (especially if 
  }						 ## Estimate has more dimensions)
}

## Need to use the Estimate array to find the
## bias, variance, and mse for each cluster
## conside looping here to get the Bias, Var, MSE for each cluster 

Bias <- apply(Estimate, 2, mean) - CM[1]  ## Modification needed here
Var <- apply(Estimate, 2, var)            ## Modification needed here
MSE <- apply((Estimate-CM[1])^2, 2, mean)  ## Modification needed here

## After you update the code, you will now want to calculate the average
## MSE of the clusters, one aveMSE for each value of lambda

## You will also want to examine the Bias, Bias^2, and Variance similarly