#Soal Penugasan modul 4

data <- c(74, 88, 93, 44, 75, 78, 92, 39, 67, 50,
          52, 43, 73, 70, 58, 78, 86, 28, 46, 61,
          63, 70, 53, 73, 80, 72, 55, 75, 65, 75,
          74, 61, 65, 38, 82, 79, 83, 59, 88, 57,
          54, 68, 65, 66, 54, 57, 62, 76, 76, 82)

# Statistik populasi
mu_pop <- mean(data)
var_pop <- var(data)

# Set parameter simulasi
set.seed(2)
ulang <- 1000
hasil <- data.frame()
# SIMULASI SAMPLING
# n = 11 sampai 20
# 1000 pengulangan per n
for (n in 11:20) {
  rata_sampel <- numeric(ulang)
  for (i in 1:ulang) {
    sampel <- sample(data, size = n, replace = TRUE)
    rata_sampel[i] <- mean(sampel)
  }
  peluang_dekat <- mean(abs(rata_sampel - mu_pop) <= 5)
  hasil <- rbind(
    hasil,
    data.frame(
      n = n,
      rata_rata_rata_sampel = mean(rata_sampel),
      sd_rata_sampel = sd(rata_sampel),
      peluang_mu_hat_dekat_mu = peluang_dekat
    )
  )
}

print("Rata-rata populasi:")
## [1] "Rata-rata populasi:"
print(mu_pop)
## [1] 66.44
print("Hasil simulasi untuk tiap n:")
## [1] "Hasil simulasi untuk tiap n:"
print(hasil)
##     n rata_rata_rata_sampel sd_rata_sampel peluang_mu_hat_dekat_mu
## 1  11              66.42991       4.346079                   0.750
## 2  12              66.37267       4.147204                   0.778
## 3  13              66.32538       4.076749                   0.781
## 4  14              66.35314       3.951366                   0.794
## 5  15              66.55093       3.803961                   0.818
## 6  16              66.48519       3.593887                   0.845
## 7  17              66.35553       3.613773                   0.838
## 8  18              66.52422       3.522843                   0.849
## 9  19              66.37584       3.339378                   0.867
## 10 20              66.54170       3.298042                   0.864
# PLOT DISTRIBUSI RATA-RATA SAMPEL
# PLOT DISTRIBUSI RATA-RATA SAMPEL
par(mfrow = c(2, 5)) # layout plot
set.seed(2)
for (n in 11:20) {
  rata_sampel <- numeric(ulang)
  for (i in 1:ulang) {
    sampel <- sample(data, size = n, replace = TRUE)
    rata_sampel[i] <- mean(sampel)
  }
  # Histogram
  hist(rata_sampel,
       main = paste("Distribusi Rata-rata Sampel (n =", n, ")"),
       xlab = "Rata-rata Sampel",
       col = "lightblue",
       freq = FALSE)
  # Kurva density dari hasil simulasi (biru)
  lines(density(rata_sampel), col = "blue", lwd = 2)
  # Kurva normal teoritis berdasarkan CLT (merah)
  curve(
    dnorm(x, mean = mean(rata_sampel), sd = sd(rata_sampel)),
    add = TRUE,
    col = "red",
    lwd = 2
  )
  # LEGEND (keterangan warna)
  legend(
    "topright",
    legend = c("Density hasil simulasi", "Kurva normal (CLT)"),
    col = c("blue", "red"),
    lwd = 2,
    bty = "n",
    cex = 0.8
  )
}