#-------------------------------------------------
#  uso di norm_mult_cond() per simulare la logica
#  inferenziale e prdittiva dietro un esperimento
#  nel quale le variabili sono distribuite normalmente
#   -> includendo ERRORI SISTEMATICI
#
#  GdA 21/2/2014
#-------------------------------------------------

# eseguire i comandi da console

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


# "esperimento" :
# la prima variabile sta per il "valore vero",
# le altre i valori misurati, distribuiti normalmente
# ad esso e "condizionatamente indipendenti"

set1 <- 2:5      # esperimenti del "set 1"
set2 <- 6:9      #   "   "      "  "set 2"
sz1  <- 0.5      # errore sistematico comune di offset del set 1
sz2  <- 1.0      # errore sistematico comune di offset del set 2
n <- max(set2)   # nr totale di osservazini sperimentali
mu <- rep(0,n)
sigma <- c(100, rep(1,n-1))  
V <- matrix(rep(sigma[1]^2, n*n), c(n,n))
V[set1, set1] <- V[set1, set1] + sz1^2    # sistematica comune set 1
V[set2, set2] <- V[set2, set2] + sz2^2    # sistematica comune set 2
diag(V)[2:n] <- diag(V)[2:n] + sigma[2:n]^2

# verifichiamo
mu                     # valori attesi
(su <- sqrt(diag(V)))  # incertezze standard  ("standard uncertainty")
V                      # matrice di covarianza
V/outer(su,su)         # matrice di covarianza

# ipotizziamo un valore della prima variabile e ricondizioniamo
mu.c <- c(2, rep(NA, n-1))
out<- norm_mult_cond(mu, V, mu.c, full=TRUE ) 
out$mu
sqrt(diag(out$V))
out$V                # le variabili > 1 sono "condizionatamente indipendenti"!

# ipotizziamo invece che la (sola) variabile 'n' assuma un preciso
# valore e ricondizioniamo
mu.c <- c(rep(NA, n-1), 2)
out<- norm_mult_cond(mu, V, mu.c, full=TRUE )
out$mu            
sqrt(diag(out$V)) 
out$V

# idem se lo facciamo sulla nr 2
mu.c <-  c(NA, 2, rep(NA, n-2))
out<- norm_mult_cond(mu, V, mu.c, full=TRUE )
out$mu
sqrt(diag(out$V))  
out$V

# assumiamo i valori del set 1 estratti random con una certa media (ad es. 5)
mu.c <- c(NA, rnorm(length(set1), 5, sigma[set1]), rep(NA, length(set2)))
out<- norm_mult_cond(mu, V, mu.c, full=TRUE )
out$mu
sqrt(diag(out$V))
out$V

# assumiamo i valori del set 2 estratti random con una certa media (ad es. 5)
mu.c <- c(NA, rep(NA, length(set1)), rnorm(length(set2), 5, sigma[set2]))
out<- norm_mult_cond(mu, V, mu.c, full=TRUE )
out$mu
sqrt(diag(out$V))
out$V


# COSA SI IMPARA?