#### Zusatzmaterial zum Normal-Normal-Modell ####

# Stichprobe erzeugen
mu <- 4
sigma2 <- 2
n <- 10
X <- rnorm(n, mean = mu, sd = sqrt(sigma2))

# Priori-Parameter
m <- 0 # Mittelwert
s <- 10 # Varianz

# Posteriori-Parameter (vgl. theoretische Ergebnisse)
mnew <-  (mean(X)*s+m*sigma2/n)/(s+sigma2/n) # Mittelwert
snew <-  (s*sigma2/n)/(s+sigma2/n) # Varianz

# Graphische Darstellung
x <- seq(-4, 8, 0.1)
plot(x, dnorm(x, mean = m, sd = sqrt(s)), type = "l", ylim = c(0,1), col = 3, xlab = "", ylab = "") # Priori
lines(x,dnorm(x, mean = mean(X), sd = sqrt(sigma2/n)), col = 1) # Likelihood
lines(x,dnorm(x, mean = mnew, sd = sqrt(snew)), col = 4) # Posteriori
legend("topleft", legend = c("Priori", "Likelihood", "Posteriori"), col = c(3,1,4), lty = 1)


# Für graphische Darstellung
x <- seq(-4, 8, 0.1)

# Verschiedene Stichprobengrößen
for(n in c(1, 5, 10, 20, 50, 100, 200, 500, 1000, 10000))
{
  # Stichprobe erzeugen
  X <- rnorm(n, mean = mu, sd = sqrt(sigma2))
  
  # Posteriori-Parameter
  mnew <-  (mean(X)*s+m*sigma2/n)/(s+sigma2/n) # Mittelwert
  snew <-  (s*sigma2/n)/(s+sigma2/n) # Varianz
  
  # Graphische Darstellung
  plot(x, dnorm(x, mean = m, sd = sqrt(s)), type = "l", ylim = c(0,1), col = 3, xlab = "", ylab = "") # Priori
  lines(x,dnorm(x, mean = mean(X), sd = sqrt(sigma2/n)), col = 1) # Likelihood
  lines(x,dnorm(x, mean = mnew, sd = sqrt(snew)), col = 4) # Posteriori
  legend("topleft", legend = c("Priori", "Likelihood", "Posteriori"), col = c(3,1,4), lty = 1)
  title(paste("n = ", n))
  
  # kurz warten
  Sys.sleep(0.5) 
}