#----------------------------------
# check di norm_mult_cond()
#  GdA 21/2/2014
#---------------------------------

# definiamo la funzione
source("norm_mult_cond.R")

# creiamo 4 variabili multivariate (con il trucco
# della trasformazione lineare, in quanto non è
# banale costruirsi matrici di covarianza per n>2)

mu0 <- c(0.4, 2, -3, 5)
sigma0 <- c(0.1, 0.2, 0.5, 2)
corr0  <- matrix(rep(0,16), c(4,4))  # tutti zeri
diag(corr0) <- rep(1,4)              # sovrascriviamo la diagonale
V0 <- outer(sigma0, sigma0) * corr0  # matrice cov diagonale

# matrice dei coefficienti 
C  <- rbind(  c(   -1,  2, -3,  1 )
             ,c(    0,  1, -1,  2 )
             ,c(    0,  3,  1, -2)
             ,c(    1, -4,  6,  0) 
           )

# valori attesi e matrice di covarianza delle combinazioni lineari
mu <- as.vector(C %*% mu0)
V <- t(C) %*% V0 %*% C
sigma <- sqrt(diag(V))
corr <- V /  outer(sigma, sigma)

eigen(V)$values # verifichiamo che la nuova V sia ok


# checks (da ESEGUIRE individualmente DA CONSOLE)
# con mu e V definite ad es in normale_multivariata.R
norm_mult_cond(mu, V, mu)
norm_mult_cond(mu, V, rep(NA,4))

norm_mult_cond(mu, V, c(NA, 21, NA, -24))
norm_mult_cond(mu, V, c(NA, 21, NA, -24), full=TRUE)

# per visualizzare la matrice di covarianza con una sola cifra decimale
round(V, 1)

# esercizio standard dottorato (-> lucidi)
mu.e <- c(2.5, 1.5)
s.e  <- c(0.5, 0.5)
rho  <- -0.8  # +0.9
V.e  <- rbind( c(s.e[1]^2, rho*s.e[1]*s.e[2])
              ,c(rho*s.e[1]*s.e[2], s.e[1]^2) )
norm_mult_cond(mu.e, V.e, c(2.0, NA))