R- Assignment 1
Thu Sep 19 2024 22:09:36 GMT+0000 (Coordinated Universal Time)
Saved by @madgakantara #r
set.seed(28092024) n <- 3 sample_norm <- rnorm(3,0,sqrt(5)) unb_estim <- var(sample_norm) max_estim <- var(sample_norm) * ((n-1)/n) mse_estim <- var(sample_norm) * ((n-1)/(n+1)) popul_var <- 5 estimates <- c(max_estim , unb_estim , mse_estim ) differences <- estimates - popul_var result <- data.frame( Estimator = c("MLE", "Unbiased", "MinMSE"), Estimate = estimates, Difference = differences ) print(result) barplot(result$Difference, names.arg = result$Estimator, main = "Difference of Estimator values fromTrue Population Variance", ylab = "Difference", col = c("red", "green", "blue"), ylim = c(min(result$Difference) - 1, max(result$Difference) + 1)) calculate_estimators <- function(n, S, mean, popul_var) { unbiased_vals <- numeric(S) mle_vals <- numeric(S) mse_vals <- numeric(S) for (i in 1:S) { sample <- rnorm(n, mean, sqrt(popul_var)) unbiased <- var(sample) mle <- unbiased * ((n - 1) / n) minmse <- unbiased * ((n - 1) / (n + 1)) unbiased_vals[i] <- unbiased mle_vals[i] <- mle mse_vals[i] <- minmse } return(list( unbiased_var = unbiased_vals, mle_var = mle_vals, mse_var = mse_vals )) } stats_calc <- function(estimate_vals, popul_var) { avg_val <- mean(estimate_vals) bias <- avg_val - popul_var variance <- var(estimate_vals) mse <- mean((estimate_vals - popul_var)^2) return(c(Average = avg_val, Bias = bias, Variance = variance, MSE = mse)) } exp_n3 = calculate_estimators(3,10000,0,sqrt(5)) exp_n20 = calculate_estimators(20,10000,0,sqrt(5)) exp_n200 = calculate_estimators(200,10000,0,sqrt(5)) results_n3 <- data.frame( Estimator = c("Unbiased", "MLE", "MinMSE"), Average = c(n3_unbiased_stats["Average"], n3_mle_stats["Average"], n3_minmse_stats["Average"]), Bias = c(n3_unbiased_stats["Bias"], n3_mle_stats["Bias"], n3_minmse_stats["Bias"]), Variance = c(n3_unbiased_stats["Variance"], n3_mle_stats["Variance"], n3_minmse_stats["Variance"]), MSE = c(n3_unbiased_stats["MSE"], n3_mle_stats["MSE"], n3_minmse_stats["MSE"]) ) results_n20 <- data.frame( Estimator = c("Unbiased", "MLE", "MinMSE"), Average = c(n20_unbiased_stats["Average"], n20_mle_stats["Average"], n20_minmse_stats["Average"]), Bias = c(n20_unbiased_stats["Bias"], n20_mle_stats["Bias"], n20_minmse_stats["Bias"]), Variance = c(n20_unbiased_stats["Variance"], n20_mle_stats["Variance"], n20_minmse_stats["Variance"]), MSE = c(n20_unbiased_stats["MSE"], n20_mle_stats["MSE"], n20_minmse_stats["MSE"]) ) results_n200 <- data.frame( Estimator = c("Unbiased", "MLE", "MinMSE"), Average = c(n200_unbiased_stats["Average"], n200_mle_stats["Average"], n200_minmse_stats["Average"]), Bias = c(n200_unbiased_stats["Bias"], n200_mle_stats["Bias"], n200_minmse_stats["Bias"]), Variance = c(n200_unbiased_stats["Variance"], n200_mle_stats["Variance"], n200_minmse_stats["Variance"]), MSE = c(n200_unbiased_stats["MSE"], n200_mle_stats["MSE"], n200_minmse_stats["MSE"]) )
Comments