simulasi monte carlo

tabel_permintaan <- data.frame(
  permintaan = c(50, 60, 70, 80, 90),
  probabilitas = c(0.1, 0.2, 0.4, 0.2, 0.1)
)
print(tabel_permintaan)

##   permintaan probabilitas
## 1         50          0.1
## 2         60          0.2
## 3         70          0.4
## 4         80          0.2
## 5         90          0.1

nilai_ekspektasi <- sum(
  tabel_permintaan$permintaan * tabel_permintaan$probabilitas
)

permintaan_simulasi <- function(n, tabel){
  kumulatif <- cumsum(tabel$probabilitas)
  u <- runif(n)
  hasil <- numeric(n)
  
  for(i in 1:n){
    for(j in 1:length(kumulatif)){
      if(u[i] <= kumulatif[j]){
        hasil[i] <- tabel$permintaan[j]
        break
      }
    }
  }
  
  return(data.frame(prediksi_permintaan = hasil))
}

set.seed(1000)
sim_1000 <- permintaan_simulasi(1000, tabel_permintaan)

cat("\n=== Frekuensi ===\n")

## 
## === Frekuensi ===

table(sim_1000$prediksi_permintaan)

## 
##  50  60  70  80  90 
## 104 188 407 205  96

mean_1000 <- mean(sim_1000$prediksi_permintaan)
cat("\nRata-rata:", mean_1000, "\n")

## 
## Rata-rata: 70.01

set.seed(5000)
mean_5000 <- mean(permintaan_simulasi(5000, tabel_permintaan)$prediksi_permintaan)

set.seed(20000)
mean_20000 <- mean(permintaan_simulasi(20000, tabel_permintaan)$prediksi_permintaan)

cat("\n=== Perbandingan ===\n")

## 
## === Perbandingan ===

cat("1000 :", mean_1000, "\n")

## 1000 : 70.01

cat("5000 :", mean_5000, "\n")

## 5000 : 69.834

cat("20000:", mean_20000, "\n")

## 20000: 69.9485

cat("\n=== Selisih dari Ekspektasi (70) ===\n")

## 
## === Selisih dari Ekspektasi (70) ===

cat("1000 :", abs(mean_1000 - nilai_ekspektasi), "\n")

## 1000 : 0.01

cat("5000 :", abs(mean_5000 - nilai_ekspektasi), "\n")

## 5000 : 0.166

cat("20000:", abs(mean_20000 - nilai_ekspektasi), "\n")

## 20000: 0.0515