Commit 33ec5cc9 authored by jbleher's avatar jbleher
Browse files

confidence bounds script added

parent 7e16332d
Loading
Loading
Loading
Loading
+485 −0

File added.

Preview size limit exceeded, changes collapsed.

+351 −0

File added.

Preview size limit exceeded, changes collapsed.

+114 −0
Original line number Diff line number Diff line
rm(list=ls())

# ------------------------------------------------------------
# Konfidenzintervalle: Frequentistische Interpretation
# ------------------------------------------------------------

set.seed(123)

mu_true <- 0         # wahrer, fester Parameter (nicht zufällig!)
sigma  <- 1          # bekannte Standardabweichung (für sauberes Z-Intervall)
n      <- 25         # Stichprobengröße
M      <- 2000       # Anzahl Wiederholungen (Stichproben)
alpha  <- 0.05
z      <- qnorm(1 - alpha/2)

# 1) Viele Stichproben ziehen und CI bauen
Xbar <- numeric(M)
lower <- numeric(M)
upper <- numeric(M)
covered <- logical(M)

for (m in 1:M) {
  XX <- rnorm(n, mean = mu_true, sd = sigma)
  Xbar[m]  <- mean(XX)
  halfwidth <- z * sigma / sqrt(n)
  lower[m] <- Xbar[m] - halfwidth
  upper[m] <- Xbar[m] + halfwidth
  covered[m] <- (lower[m] <= mu_true) && (mu_true <= upper[m])
}

# 2) Empirische Überdeckungswahrscheinlichkeit
cat("Empirische Überdeckung (sollte ~ 0.95 sein): ",
    mean(covered), "\n")

# 3) Plot: viele Intervalle, rot wenn verfehlt, schwarz wenn getroffen
K <- 100  # wie viele Intervalle im Plot gezeigt werden
idx <- sample(1:M, K)

par(mar = c(5, 4, 6, 2) + 0.1)  # bottom, left, top, right
plot(NULL, xlim = range(c(lower[idx], upper[idx])), ylim = c(1, K),
     xlab = expression(mu), ylab = "Intervall-Index",
     main = "95%-Konfidenzintervalle: Viele Stichproben, fixer Parameter")
abline(v = mu_true, lwd = 2, lty = 2)  # wahrer Parameter (fix)

for (i in 1:K) {
  j <- idx[i]
  col <- if (covered[j]) "black" else "red"
  segments(lower[j], i, upper[j], i, col = col, lwd = 2)
}

legend("top",
       inset = c(0, -0.07),
       legend = c("enthält mu (Treffer)", "enthält mu nicht (Fehler)", "wahrer mu"),
       col = c("black", "red", "black"),
       lty = c(1, 1, 2),
       lwd = c(2, 2, 2),
       bty = "n",
       horiz = TRUE,
       xpd = NA)

# ------------------------------------------------------------
# 4) EIN beobachtetes Konfidenzintervall – und was wir sinnvoll
#    darüber sagen können (frequentistische Sicht)
# ------------------------------------------------------------

# Eine konkrete (einmalige) Stichprobe
x_obs <- rnorm(n, mean = mu_true, sd = sigma)
xbar_obs <- mean(x_obs)
halfwidth_obs <- z * sigma / sqrt(n)
ci_obs <- c(xbar_obs - halfwidth_obs, xbar_obs + halfwidth_obs)

cat("\nEin konkret beobachtetes 95%-Konfidenzintervall:\n")
cat("[", round(ci_obs[1], 3), ", ", round(ci_obs[2], 3), "]\n", sep = "")

# Wichtige Klarstellung:
cat("\nWICHTIG:\n")
cat("Der wahre Parameter mu ist unbekannt.\n")
cat("Wir wissen NICHT, ob dieses konkrete Intervall mu enthält.\n")
cat("Die Zahl 95% ist KEINE Wahrscheinlichkeit für dieses Intervall.\n")

# ------------------------------------------------------------
# 4a) Was bedeutet die 95%-Aussage dann?
# ------------------------------------------------------------

cat("\nInterpretation der 95%:\n")
cat("Das verwendete Intervall-Verfahren hat eine Überdeckungswahrscheinlichkeit von 95%.\n")
cat("Das heißt:\n")
cat("- Würden wir dieses Experiment sehr oft wiederholen,\n")
cat("- dann würden ca. 95% der resultierenden Intervalle den wahren Parameter enthalten.\n")

# ------------------------------------------------------------
# 4b) Gedankliches Experiment: Wiederholung des Verfahrens
#      (mu bleibt fix, Stichprobe variiert)
# ------------------------------------------------------------

M_long <- 10000
covered_long <- logical(M_long)

for (m in 1:M_long) {
  x <- rnorm(n, mean = mu_true, sd = sigma)
  xbar <- mean(x)
  ci <- c(xbar - z * sigma / sqrt(n),
          xbar + z * sigma / sqrt(n))
  covered_long[m] <- (ci[1] <= mu_true) && (mu_true <= ci[2])
}

cat("\nLangfristige Überdeckungsrate des Verfahrens:\n")
cat("Anteil Intervalle mit mu enthalten = ",
    round(mean(covered_long), 3), "\n", sep = "")


cat("\nKernaussage:\n")
cat("Ein 95%-Konfidenzintervall ist eine Aussage über die ZUVERLÄSSIGKEIT\n")
cat("des Verfahrens, nicht über die Wahrscheinlichkeit eines einzelnen Intervalls.\n")