% This code implements the Rosenblatt Kernel density estimator
% data generation from a mixture of gaussian distributions or a single
% gaussian

% Sample size
n=2500;

% Random number generation (gaussian)

X=randn(n,1);
f=(1/sqrt(2*pi)).*exp(-0.5*X.^2);

% Random number generation (mixed gaussian)
%[x, f] = mixnormtrimodal(n);


% Bandwidth selection by Cross Validation

% initial value
c=.2;
hpi = n^(-1/5)*std(X)*c;

% The next line calls the function fminsearch with function handle
% ugcv(x,h)

[hmin,value,flag,output]=fminsearch(@(h)ugcv(X,h),hpi);

%The estimator - calls the function urosenblatt
fhat = urosenblatt(X,X,hmin,5);
 
[hmin hpi]
% Prepare to plot

xf = [X f fhat];
xf=sortrows(xf,1);

figure(1)
plot(xf(:,1),xf(:,2),xf(:,1),xf(:,3)), title('Density estimation')
box off
grid on
xlabel('Values of x')
ylabel('Density, Estimated density')

