#------------------------#
#                        #
#   Loesungen Blatt 8    #
#                        #
#------------------------#



# Aufgabe 1)


# 1a)

N = 10000

mu1 = 5
mu2 = -2
mu3 = -2
mu4 = 0

sigma1 = sqrt(3)
sigma2 = 1 
sigma3 = 1 
sigma4 = 2

phi1 = rnorm(N, mean=mu1, sd=sigma1 )
phi2 = rnorm(N, mean=mu2, sd=sigma2 )
phi3 = rnorm(N, mean=mu3, sd=sigma3 )
phi4 = rnorm(N, mean=mu4, sd=sigma4 )

phi = phi1 + phi2 + phi3 + phi4
phi                                     # ok, ist Vektor mit Laenge N=10000


# 1b)

hist(phi, prob=TRUE, breaks=40 )


# 1c)

mu = mu1 + mu2 + mu3 + mu4
sigma = sqrt( sigma1^2 + sigma2^2 + sigma3^2 + sigma4^2 )
mu
sigma

curve( dnorm(x, mean=mu, sd=sigma), col="red", add=TRUE)        # das passt..



# Aufgabe 2)


# 2a)

curve( pnorm(x), from=-5, to=+5 )
curve( pnorm(-x), add=TRUE, col="green")
curve( pnorm(x)+pnorm(-x), add=TRUE, col="red")


# 2b)  ist im Loesung8.pdf




# Aufgabe 3) 


# 3a)

kk = 1:6

prob = pnorm(kk) - pnorm(-kk)
prob

options(digits = 10)          # mehr Stellen, die letzte W'keit ist nicht wirklich = 1
prob


# 3b) 

# sind dieselben Zahlen wie in 3a, siehe auch Loesung8.pdf


# 3c)

# Wir wollen haben:       N * Prob[ |y| <= k ] = N - 1 
#
#               <=>        Prob[ |y| <= k ] = 1 - 1/N 
#
# Es ist               Prob[ |y| <= k ] = pnorm(k) - pnorm(-k) 
#
# und wegen (Aufg2)           pnorm(k) + pnorm(-k) = 1   
#
# ist                         pnorm(-k) = 1 - pnorm(k) 
#
# Also: 
#             Prob[ |y| <= k ] = pnorm(k) - pnorm(-k) = 2 * pnorm(k) - 1 
#
# und das soll gleich    1 - 1/N    sein:
#
#                         2 * pnorm(k) - 1 = 1 - 1/N    
#
#                   <=>   2 * pnorm(k)     = 2 - 1/N    
#
#                   <=>       pnorm(k)     = 1 - 1/(2N) 
#
# Um das k zu finden, brauchen wir jetzt die Umkehrfunktion von dem pnorm, 
# das war gerade das qnorm, die Quantil-Funktion. Also bekommen wir 
#
#                           k = qnorm( 1 - 1/(2N) )
#
# Geben wir noch die numerischen Werte an:

N = 10^6

k = qnorm( 1 - 1/(2*N) )
k                            # etwa 4.89

N = 10^9

k = qnorm( 1 - 1/(2*N) )
k                            # etwa 6.11


# Schauen wir uns das noch konkret mit 1 Mio Zufallszahlen an:

y  = rnorm(10^6)             # 1 Mio standard-normalverteilte Zufallszahlen 
absy = abs(y)                # wir nehmen den Betrag
sy = sort(absy)              # jetzt der Groesse nach geordnet, von klein nach gross

tail(sy, 10)                 # die groessten 10 Zahlen, die letzten ein, zwei oder drei dann so 
                             # typischerweise immer ein bischen ueber oder unter 4.89, das passt 
                             # so im wesentlichen..