# SOL Lab 7 TSdD 2024

# Funzioni di perdita e criteri di ottimalità


# Es.1
# Confronto decisioni con VA, mM, SC

# Grafici
# (a)
W.1.fun=function(x){x^2}
W.2.fun=function(x){1-x^2}
W.3.fun=function(x){0.5*x}

# (b)
curve(W.1.fun,xlim=c(0,1),ylab="perdite", xlab="omega",lty=3,lwd=2)
curve(W.2.fun,add=T,lty=3)
curve(W.3.fun,add=T,lty=4)

# (c)
lam=1/4
abline(h=lam)
title("Funzioni di Perdita")

(d)
# Simulazione da p(w)
M=5000
w.MC=rbeta(M,1,1)

# lam=1/8 # lambda del K.sc (soglia critica)

# (e)
W.1=w.MC^2
W.2=1-w.MC
W.3=0.5*w.MC

# (f)
# Criteri di ottimalità
K.va.1=mean(W.1)
K.va.2=mean(W.2)
K.va.3=mean(W.3)
K.va.1
K.va.2
K.va.2


K.sc.1=mean(W.1>lam)
K.sc.2=mean(W.2>lam)
K.sc.3=mean(W.3>lam)
K.sc.1
K.sc.2
K.sc.3

# (g)
lam=1/2
K.sc.1=mean(W.1>lam)
K.sc.2=mean(W.2>lam)
K.sc.3=mean(W.3>lam)
K.sc.1
K.sc.2
K.sc.3

# (h)
w.MC=rbeta(M,1,4)    # e ripetere




# Es. 2

# (a)
delta=1/2
W.quad.fun=function(w){(w-delta)^2}

W.quad.gen.fun=function(w){
  ((w-delta)^2)/(w*(1-w))
}

W.abs.fun=function(w){
  abs(w-delta)
}

# (b)

curve(W.quad.fun(x),xlim=c(0,1),lty=3,lwd=2,ylab="perdite",xlab="omega",ylim=c(0,1))
curve(W.quad.gen.fun(x),lty=2, add=T,lwd=2)
curve(W.abs.fun(x),lty=4, add=T,lwd=2)

# (c)
# guardare i grafici e usare le funzioni

aa=seq(0,1,0.1)

max(W.quad.fun(aa))
W.quad.fun(0)

max(W.quad.gen.fun(aa))
W.quad.gen.fun(0)

max(W.abs.fun(aa))
W.abs.fun(0)



# (d)
# Simulazioni

M=5000
w.MC=rbeta(M,2,2)

# (e)
mean(w.MC)


# (f)
W.quad.MC=(w.MC-delta)^2
W.quad.gen.MC=((w.MC-delta)^2)/(w.MC*(1-w.MC))
W.abs.MC=abs(w.MC-delta)

# Criterio VA
K.va.quad=mean(W.quad.MC)
K.va.quad.gen=mean(W.quad.gen.MC)
K.va.abs=mean(W.abs.MC)
K.va.quad
K.va.quad.gen
K.va.abs

# (g)
# Criterio SC (soglia critica)
lam=1/8 
K.sc.quad=mean(W.quad.MC>lam)
K.sc.quad.gen=mean(W.quad.gen.MC>lam)
K.sc.abs=mean(W.abs.MC>lam)
K.sc.quad
K.sc.quad.gen
K.sc.abs


# Es. 3

M.omega=1000
omega=runif(M.omega)
omega=matrix(omega,1,length(omega))
# omega

# (a)
dd=seq(0,1,.1)
dd=matrix(dd,length(dd),1)
dd
N.d=length(dd)
N.d

# (b)
matrix.d.o=matrix(0,N.d,M.omega)

for(i in 1:N.d){
  for(j in 1:M.omega){
    matrix.d.o[i,j]=(dd[i]-omega[j])^2
  }}

dim(matrix.d.o)

# (c)
E.W=apply(matrix.d.o,1,mean)
length(E.W)

plot(dd,E.W,type="l")

# (e)
dd[E.W==min(E.W)]