function M=nmdmise(h,n,w,sigma,mu)

% NMDMISE     MISE(h) for normal mixture densities.
%
%             Example: To reproduce Fig. 2.5 of Wand and Jones use:
%                      
%             for i=0:100,
%               h(i+1)=0.5+(i-40)/100; 
%               M(i+1)=nmdmise(h(i+1),100,[0.75 0.25]',[1 1/3]',[0 3/2]'); 
%             end; 
%             plot(h,M);
%             [Mval,ind]=min(M);
%             h_min=h(ind)
%
%             NMDMISE(H,N,W,SIGMA,MU)
%
%             Christian C. Beardah 1995 

k=length(w);

M=1/(2*sqrt(pi)*n*h);

for i=1:k,
  for j=1:k,
    omega0(i,j)=npdf(mu(i)-mu(j),0,sqrt(sigma(i)^2+sigma(j)^2));
    omega1(i,j)=npdf(mu(i)-mu(j),0,sqrt(h^2+sigma(i)^2+sigma(j)^2));
    omega2(i,j)=npdf(mu(i)-mu(j),0,sqrt(2*h^2+sigma(i)^2+sigma(j)^2));    
  end;
end;

M=M+w'*((1-1/n)*omega2-2*omega1+omega0)*w;