## ESTIMADOR NAIVE

uni<-function(x){apply(x,2,function(x) ifelse(abs(x)<1,1/2,0))}

est_nucleo<-function(t,x,bandwith,nucleo){
apply(nucleo(outer(t,x,"-")/bandwith),1,"mean")/bandwith}


## DATOS
n<-100
M=100
bandwith<-seq(0,8,by=0.05)
t<-seq(-3,3,by=0.05)
muestra<-rnorm(n)
F<-dnorm(t)

## MISE es calculado 100 veces para cada h. 
MISE<-matrix(0,nrow=M,ncol=length(bandwith))

for(m in 1:M){
MISE[m,]<-apply(as.matrix(bandwith),1,function(h){sum((est_nucleo(t,muestra,h,uni)-F)^2)})
}

Sol<-apply(MISE,2,"mean")
plot(Sol~bandwith, type="l")
h_opt<-bandwith[which(Sol==min(Sol,na.rm=T))]
