# LDS - Prova finale 1 - dicembre 2024 - Compito A

# Es. 1
rm(list = ls())
al=2
be=3
n=10
sn=4
M=50000
al.p=al+n
be.p=be+sn
theta.mc=rgamma(M,al.p,rate=be.p)
psi=1/theta.mc^2
mean(psi)
# 0.445


# Es. 2
rm(list = ls())
al=3
be=4
n=10
sn=3
al.p=al+sn
be.p=be+n
gg=1-0.85
L=qgamma(gg/2,al.p,be.p)
U=qgamma(1-gg/2,al.p,be.p)
c(exp(-U),exp(-L))
# (0,496,0.812)


# Es. 3
rm(list = ls())
delta=2
M=50000
lambda=2
omega.mc=rnorm(M,mean=4,sd=sqrt(2))
W=abs(delta-omega.mc)
mean(W)
mean(W>lambda)
# 0.506


# Es. 4
rm(list = ls())
# rischio R coincide con varianza dello stimatore, pari a theta(1-theta)/n
n=5
al=5
be=3
M=50000
theta.mc=rbeta(M,al,be)
R=theta.mc*(1-theta.mc)/n
# r=E(R)
r=mean(R)
r
# 0.041


# Es. 5
rm(list = ls())
mu=5
sig=sqrt(2)
n=10
M=50000
x.matr=rnorm(n*M,mean=mu,sd=sig)
x.matr=matrix(x.matr,M,n)
x.matr3=x.matr^3
m3=apply(x.matr3,1,mean)
Ex.m3=mean(m3)
Ex.m3
# 155
# Se fatto con m2
x.matr=rnorm(n*M,mean=mu,sd=sig)
x.matr=matrix(x.matr,M,n)
x.matr2=x.matr^2
m2=apply(x.matr2,1,mean)
Ex.m2=mean(m2)
Ex.m2
# 935.5


# Es. 6
rm(list = ls())
M=50000
theta.vero=3
n=15
x.matr=rgamma(M*n,1,rate=theta.vero)
# x.matr=rexp(M*n,theta.vero)
x.matr=matrix(x.matr,M,n)
sn.mc=apply(x.matr,1,sum)
L=n/sn.mc - n/(sn.mc*sqrt(n))
U=n/sn.mc + n/(sn.mc*sqrt(n))
# copertura
mean(L<theta.vero&U>theta.vero)
# 0.68


# Es. 7
rm(list = ls())
M=50000
th.vero=3
n=9
k=1
# Per il calcolo di alpha (simulare sotto H0)
x.matr=rnorm(n*M,mean=0,sd=sqrt(th.vero))
x.matr=matrix(x.matr,M,n)
medie=apply(x.matr,1,mean)
for(i in 1:M){
  x.matr[i,]=abs(x.matr[i,]-medie[i])
}
A=apply(x.matr,1,mean)
mean(A<k)
# 0.19



# Es. 8
rm(list = ls())
M=50000
th.0=2
th.1=1
th.vero=1.5
x.mc=rgamma(M,1,th.vero)
x.mc=matrix(x.mc,M,1)
NUM.0=dgamma(x.mc,1,rate=th.0)
DEN.1=dgamma(x.mc,1,rate=th.1)
lambda=NUM.0/DEN.1
E=mean(lambda)
E
# 1.2