mu <- c(0, 2, 5, 10, 1.5)      # Nilai Rata-rata (mu)
sigma <- c(1, 1, 5, 3, 0.1)    # Nilai Simpangan Baku (sigma)

x <- seq(-5, 20, length.out = 500)

plot(x, dnorm(x, mean = mu[1], sd = sigma[1]), type = "n", 
     main = "Simulasi Kurva Distribusi Normal Tertumpuk",
     xlab = "X", ylab = "Densitas", ylim = c(0, 4), xlim = c(-4, 18))

warna <- c("blue", "red", "green", "purple", "darkorange")

for (i in 1:length(mu)) {
  y <- dnorm(x, mean = mu[i], sd = sigma[i])
  lines(x, y, col = warna[i], lwd = 2)
  polygon(c(x[1], x, x[length(x)]), c(0, y, 0), 
          col = adjustcolor(warna[i], alpha.f = 0.2), 
          border = warna[i])
}

x <- seq(-4, 4, length.out = 500)
df_values <- c(1, 3, 10, 30)

plot(x, dt(x, df = df_values[1]), type = "n", 
     main = "Simulasi Kurva Distribusi t Tertumpuk",
     xlab = "Nilai t", ylab = "Densitas", ylim = c(0, 0.4))

warna <- c("blue", "red", "green", "purple")

for (i in 1:length(df_values)) {
  y <- dt(x, df = df_values[i])
  lines(x, y, col = warna[i], lwd = 2)
  polygon(c(x[1], x, x[length(x)]), c(0, y, 0), 
          col = adjustcolor(warna[i], alpha.f = 0.2), 
          border = warna[i])
}

legend("topright", 
       legend = paste("t-dist (df =", df_values, ")"), 
       fill = adjustcolor(warna, alpha.f = 0.2), 
       border = warna)

#distribusi f

x_f <- seq(0, 5, length.out = 500)
df1_vals <- c(2, 5, 10)
df2_vals <- c(5, 2, 10)
warna_f <- c("blue", "red", "green")

plot(x_f, df(x_f, df1 = df1_vals[1], df2 = df2_vals[1]), type = "n", 
     main = "Simulasi Kurva Distribusi F",
     xlab = "Nilai F", ylab = "Densitas", ylim = c(0, 1))

for (i in 1:length(df1_vals)) {
  y_f <- df(x_f, df1 = df1_vals[i], df2 = df2_vals[i])
  lines(x_f, y_f, col = warna_f[i], lwd = 2)
  polygon(c(x_f[1], x_f, x_f[length(x_f)]), c(0, y_f, 0), 
          col = adjustcolor(warna_f[i], alpha.f = 0.2), border = warna_f[i])
}

legend("topright", legend = paste("F-dist (df1=", df1_vals, ", df2=", df2_vals, ")"), 
       fill = adjustcolor(warna_f, alpha.f = 0.2), border = warna_f, cex = 0.8)

#distribusi chi-square

x_chi <- seq(0, 20, length.out = 500)
df_chi <- c(2, 4, 8)
warna_chi <- c("purple", "orange", "brown")

plot(x_chi, dchisq(x_chi, df = df_chi[1]), type = "n", 
     main = "Simulasi Kurva Distribusi Chi-Square",
     xlab = "Nilai Chi-Square", ylab = "Densitas", ylim = c(0, 0.5))

for (i in 1:length(df_chi)) {
  y_chi <- dchisq(x_chi, df = df_chi[i])
  lines(x_chi, y_chi, col = warna_chi[i], lwd = 2)
  polygon(c(x_chi[1], x_chi, x_chi[length(x_chi)]), c(0, y_chi, 0), 
          col = adjustcolor(warna_chi[i], alpha.f = 0.2), border = warna_chi[i])
}

legend("topright", legend = paste("Chi-Sq (df =", df_chi, ")"), 
       fill = adjustcolor(warna_chi, alpha.f = 0.2), border = warna_chi, cex = 0.8)

#distribusi gamma

x_gamma <- seq(0, 15, length.out = 500)
shape_vals <- c(1, 2, 3)
rate_vals <- c(1, 1, 1)
warna_gamma <- c("blue", "red", "green")

plot(x_gamma, dgamma(x_gamma, shape = shape_vals[1], rate = rate_vals[1]), type = "n", 
     main = "Simulasi Kurva Distribusi Gamma",
     xlab = "x", ylab = "Densitas", ylim = c(0, 1))

for (i in 1:length(shape_vals)) {
  y_gamma <- dgamma(x_gamma, shape = shape_vals[i], rate = rate_vals[i])
  lines(x_gamma, y_gamma, col = warna_gamma[i], lwd = 2)
  polygon(c(x_gamma[1], x_gamma, x_gamma[length(x_gamma)]), c(0, y_gamma, 0), 
          col = adjustcolor(warna_gamma[i], alpha.f = 0.2), border = warna_gamma[i])
}

legend("topright", legend = paste("Gamma (shape=", shape_vals, ")"), 
       fill = adjustcolor(warna_gamma, alpha.f = 0.2), border = warna_gamma, cex = 0.8)

#distribusi weibul

x_weibull <- seq(0, 3, length.out = 500)
shape_weib <- c(1, 2, 5)
scale_weib <- c(1, 1, 1)
warna_weib <- c("purple", "orange", "brown")

plot(x_weibull, dweibull(x_weibull, shape = shape_weib[1], scale = scale_weib[1]), type = "n", 
     main = "Simulasi Kurva Distribusi Weibull",
     xlab = "x", ylab = "Densitas", ylim = c(0, 2.5))

for (i in 1:length(shape_weib)) {
  y_weib <- dweibull(x_weibull, shape = shape_weib[i], scale = scale_weib[i])
  lines(x_weibull, y_weib, col = warna_weib[i], lwd = 2)
  polygon(c(x_weibull[1], x_weibull, x_weibull[length(x_weibull)]), c(0, y_weib, 0), 
          col = adjustcolor(warna_weib[i], alpha.f = 0.2), border = warna_weib[i])
}

legend("topright", legend = paste("Weibull (shape=", shape_weib, ")"), 
       fill = adjustcolor(warna_weib, alpha.f = 0.2), border = warna_weib, cex = 0.8)

#distribusi eksponensial

x_exp <- seq(0, 5, length.out = 500)
rate_exp <- c(0.5, 1, 2)
warna_exp <- c("darkblue", "darkred", "darkgreen")

plot(x_exp, dexp(x_exp, rate = rate_exp[1]), type = "n", 
     main = "Simulasi Kurva Distribusi Eksponensial",
     xlab = "x", ylab = "Densitas", ylim = c(0, 2))

for (i in 1:length(rate_exp)) {
  y_exp <- dexp(x_exp, rate = rate_exp[i])
  lines(x_exp, y_exp, col = warna_exp[i], lwd = 2)
  polygon(c(x_exp[1], x_exp, x_exp[length(x_exp)]), c(0, y_exp, 0), 
          col = adjustcolor(warna_exp[i], alpha.f = 0.2), border = warna_exp[i])
}

legend("topright", legend = paste("Exp (rate=", rate_exp, ")"), 
       fill = adjustcolor(warna_exp, alpha.f = 0.2), border = warna_exp, cex = 0.8)

#distribusi pareto

x_pareto <- seq(1, 10, length.out = 500)
xm <- 1
alpha_vals <- c(1, 2, 3)
warna_par <- c("blue", "red", "green")

# Fungsi PDF Pareto manual: alpha * xm^alpha / x^(alpha + 1)
dpareto_manual <- function(x, xm, alpha) {
  ifelse(x >= xm, (alpha * xm^alpha) / (x^(alpha + 1)), 0)
}

plot(x_pareto, dpareto_manual(x_pareto, xm, alpha_vals[1]), type = "n", 
     main = "Simulasi Kurva Distribusi Pareto",
     xlab = "x", ylab = "Densitas", ylim = c(0, 3))

for (i in 1:length(alpha_vals)) {
  y_par <- dpareto_manual(x_pareto, xm, alpha_vals[i])
  lines(x_pareto, y_par, col = warna_par[i], lwd = 2)
  polygon(c(x_pareto[1], x_pareto, x_pareto[length(x_pareto)]), c(0, y_par, 0), 
          col = adjustcolor(warna_par[i], alpha.f = 0.2), border = warna_par[i])
}

legend("topright", legend = paste("Pareto (alpha=", alpha_vals, ")"), 
       fill = adjustcolor(warna_par, alpha.f = 0.2), border = warna_par, cex = 0.8)

#distribusi uniform

x_unif <- seq(-1, 6, length.out = 500)
min_vals <- c(0, 1, 0)
max_vals <- c(3, 4, 5)
warna_unif <- c("purple", "orange", "darkcyan")

plot(x_unif, dunif(x_unif, min = min_vals[1], max = max_vals[1]), type = "n", 
     main = "Simulasi Kurva Distribusi Uniform",
     xlab = "x", ylab = "Densitas", ylim = c(0, 0.6))

for (i in 1:length(min_vals)) {
  y_unif <- dunif(x_unif, min = min_vals[i], max = max_vals[i])
  lines(x_unif, y_unif, col = warna_unif[i], lwd = 2)
  polygon(c(x_unif[1], x_unif, x_unif[length(x_unif)]), c(0, y_unif, 0), 
          col = adjustcolor(warna_unif[i], alpha.f = 0.2), border = warna_unif[i])
}

legend("topright", legend = paste("Uniform (min=", min_vals, ", max=", max_vals, ")"), 
       fill = adjustcolor(warna_unif, alpha.f = 0.2), border = warna_unif, cex = 0.8)

x_beta <- seq(0, 1, length.out = 500)
shape1_vals <- c(1, 2, 2, 5)
shape2_vals <- c(1, 2, 5, 2)
warna_beta <- c("blue", "red", "green", "purple")

plot(x_beta, dbeta(x_beta, shape1 = shape1_vals[1], shape2 = shape2_vals[1]), type = "n", 
     main = "Simulasi Kurva Distribusi Beta",
     xlab = "x", ylab = "Densitas", ylim = c(0, 3))

for (i in 1:length(shape1_vals)) {
  y_beta <- dbeta(x_beta, shape1 = shape1_vals[i], shape2 = shape2_vals[i])
  lines(x_beta, y_beta, col = warna_beta[i], lwd = 2)
  polygon(c(x_beta[1], x_beta, x_beta[length(x_beta)]), c(0, y_beta, 0), 
          col = adjustcolor(warna_beta[i], alpha.f = 0.2), border = warna_beta[i])
}

legend("top", legend = paste("Beta (a=", shape1_vals, ", b=", shape2_vals, ")"), 
       fill = adjustcolor(warna_beta, alpha.f = 0.2), border = warna_beta, cex = 0.7)