#Distribuciķ Normal
#Evoluciķ de l'amplada de l'interval de confianįa per a la Proporiciķ de la poblaciķ
#a mesura que incrementa n i per diferents alfas

#Fixem la llavor per a que l'experiment sempre surti el mateix

set.seed(1000)

#Fiquem les dades
p<-0.4
n<-1
rep<-500
nvars<-seq(30,100,1)
mu<-p

#Creem els vectors necessaris de la Mitjana Mostral
e_mu0.hat<-vector()
eII_mu0.hat<-vector()
eIS_mu0.hat<-vector()
e2_mu0.hat<-vector()
eII2_mu0.hat<-vector()
eIS2_mu0.hat<-vector()

e_mu1.hat<-vector()
eII_mu1.hat<-vector()
eIS_mu1.hat<-vector()
e2_mu1.hat<-vector()
eII2_mu1.hat<-vector()
eIS2_mu1.hat<-vector()

e_mu2.hat<-vector()
eII_mu2.hat<-vector()
eIS_mu2.hat<-vector()
e2_mu2.hat<-vector()
eII2_mu2.hat<-vector()
eIS2_mu2.hat<-vector()

for (i in 1:length(nvars)){
	matriu.x<-matrix(rbinom(rep*nvars[i],n,p),rep)
	sigma<-sqrt(p*(1-p))
#Interval alfa
	alfa=0.1
	mu0.hat<-rowMeans(matriu.x)
	II_mu0.hat<-mu0.hat+qnorm(alfa/2)*sigma/sqrt(nvars[i])
	IS_mu0.hat<-mu0.hat+qnorm(1-alfa/2)*sigma/sqrt(nvars[i])
	e_mu0.hat[i]<-mean(mu0.hat)
	eII_mu0.hat[i]<-mean(II_mu0.hat)
	eIS_mu0.hat[i]<-mean(IS_mu0.hat)
	e2_mu0.hat<-rbind(e2_mu0.hat,e_mu0.hat[i])
	eII2_mu0.hat<-rbind(eII2_mu0.hat,eII_mu0.hat[i])
	eIS2_mu0.hat<-rbind(eIS2_mu0.hat,eIS_mu0.hat[i])
	alfa=0.05
	mu1.hat<-rowMeans(matriu.x)
	II_mu1.hat<-mu1.hat+qnorm(alfa/2)*sigma/sqrt(nvars[i])
	IS_mu1.hat<-mu1.hat+qnorm(1-alfa/2)*sigma/sqrt(nvars[i])
	e_mu1.hat[i]<-mean(mu1.hat)
	eII_mu1.hat[i]<-mean(II_mu1.hat)
	eIS_mu1.hat[i]<-mean(IS_mu1.hat)
	e2_mu1.hat<-rbind(e2_mu1.hat,e_mu1.hat[i])
	eII2_mu1.hat<-rbind(eII2_mu1.hat,eII_mu1.hat[i])
	eIS2_mu1.hat<-rbind(eIS2_mu1.hat,eIS_mu1.hat[i])
	alfa=0.01
	mu2.hat<-rowMeans(matriu.x)
	II_mu2.hat<-mu2.hat+qnorm(alfa/2)*sigma/sqrt(nvars[i])
	IS_mu2.hat<-mu2.hat+qnorm(1-alfa/2)*sigma/sqrt(nvars[i])
	e_mu2.hat[i]<-mean(mu2.hat)
	eII_mu2.hat[i]<-mean(II_mu2.hat)
	eIS_mu2.hat[i]<-mean(IS_mu2.hat)
	e2_mu2.hat<-rbind(e2_mu2.hat,e_mu2.hat[i])
	eII2_mu2.hat<-rbind(eII2_mu2.hat,eII_mu2.hat[i])
	eIS2_mu2.hat<-rbind(eIS2_mu2.hat,eIS_mu2.hat[i])
}
res0<-cbind(e2_mu0.hat,eII2_mu0.hat,eIS2_mu0.hat)
res1<-cbind(e2_mu1.hat,eII2_mu1.hat,eIS2_mu1.hat)
res2<-cbind(e2_mu2.hat,eII2_mu2.hat,eIS2_mu2.hat)

nomscol<-c("Esperanįa", "I. Inferior", "I. Superior")
colnames(res0)<-nomscol
colnames(res1)<-nomscol
colnames(res2)<-nomscol

res0; res1; res2; 

#Anālisi grāfica.
#Generem una matriu grāfica d'una fila i dues columnes
par(mfrow = c(1,3), oma=c(1, 0, 4, 0))

#Grāfic a la posiciķ 1,1
y<-seq(mu-0.75*sigma, mu+0.75*sigma, along=nvars)
x<-nvars
plot(x,y, main='Alfa=0.1',xlab = "Mida mostral", ylab = "", type="n")
abline(h=mu)
lines(x,e2_mu0.hat, lty=1, lwd=3, col='red')
lines(x,eII2_mu0.hat, lty=1, lwd=3, col='blue')
lines(x,eIS2_mu0.hat, lty=1, lwd=3, col='blue')


#Grāfic a la posiciķ 1,2
y<-seq(mu-0.75*sigma, mu+0.75*sigma, along=nvars)
x<-nvars
plot(x,y, main='Alfa=0.05',xlab = "Mida mostral", ylab = "", type="n")
abline(h=mu)
lines(x,e2_mu1.hat, lty=1, lwd=3, col='red')
lines(x,eII2_mu1.hat, lty=1, lwd=3, col='blue')
lines(x,eIS2_mu1.hat, lty=1, lwd=3, col='blue')

#Grāfic a la posiciķ 1,3
y<-seq(mu-0.75*sigma, mu+0.75*sigma, along=nvars)
x<-nvars
plot(x,y, main='Alfa=0.01',xlab = "Mida mostral", ylab = "", type="n")
abline(h=mu)
lines(x,e2_mu2.hat, lty=1, lwd=3, col='red')
lines(x,eII2_mu2.hat, lty=1, lwd=3, col='blue')
lines(x,eIS2_mu2.hat, lty=1, lwd=3, col='blue')
mtext(side=3, line=0, cex=1.5, outer=T,"Interval de confianįa per a la Proporciķ Poblacional. Variāncia coneguda.")
mtext(side=1, line=-1, cex=0.75, outer=T, adj=1, "Script creat per Jordi Lķpez-Tamayo  ")