#-------------------------------------------------------------------
# prova di importance sampling:
#
#  vogliamo stimare il valore medio di exp(-abs(x)/tau)
#  con x distribuito secondo una gaussiana di parametri mu.p e sigma.p
#
# GdA 15 marzo 2012
##-------------------------------------------------------------------- 

n <- 100000 # numero di eventi
s = 1000    # step di controllo intermedio

# parametro di "f(x)"
tau=0.2

# parametri gaussiana "p(x)"
mu.p <- 3
sigma.p <- 1

# parametri gaussiana "q(x)"
mu.q <- 0
sigma.q <- 1

# campionamento normale
x <- rnorm(n, mu.p, sigma.p)
fx <- exp(-abs(x)/tau)
ni <- seq(s, n, s)
Ifx = rep(0, length(ni))
for (i in 1:length(ni)) {
  Ifx[i] = sum(fx[1:ni[i]])/ni[i]
}

# importance sampling
x <- rnorm(n, mu.q, sigma.q)
fx <- exp(-abs(x)/tau)*dnorm(x, mu.p, sigma.p)/dnorm(x, mu.q, sigma.q)
Ifx.IS = rep(0, length(ni))
for (i in 1:length(ni)) {
  Ifx.IS[i] = sum(fx[1:ni[i]])/ni[i]
}

plot(ni,Ifx, ylim=c(min(Ifx,Ifx.IS),max(Ifx,Ifx.IS)), col='blue', ty='l')
points(ni,Ifx.IS, ty='l', col='red')