Usage: 
{f,band}=EBBSmain(x,y,xgrid,hgrid,OrderDer,p,type,varest{,msespan{,nterms{,bandspan{,Kernel{,J2{,J1{,RidgeCoef}}}}}}})

Input: 
 x  n x k matrix, the independent variables.

 y  n x 1 vector, the dependent variable.

 xgrid  m x k matrix, grid where you
estimate the dependent variable.

 hgrid  l x 1 vector, grid of bandwidths. For each point of xgrid, this
algorithm determines the best bandwidth in hgrid. Take care that
the first J1 values and the last J2 values are excluded.
Moreover, as the data are standardized, the grid is always a
vector even if the data are multivariate.

 OrderDer  scalar, order of the derivative you want to estimate. This order
depends on the dimension of the independent variables.
 OrderDer = 0 corresponds always to the function itself
 1 <= OrderDer <= k corresponds to the first derivatives
 k < OrderDer <= 2*k corresponds to the non mixed
second derivatives
 OrderDer > 2*k corresponds to the mixed derivatives

 p  scalar, degree of the polynomials (1<=p<=2). Take care that if you want to
estimate second derivatives, then p must be equal to 2.

 type  scalar, type of variance estimation:
 1 means that "varest" contains a variance estimation of the
response for each target point in xgrid.
 2 means that the algorithm will use "coef=(sqrt(var(y))/mean(y))^2"
to estimate the variance part in the MSE. Please refer to Ruppert's article.

 varest  if type==1, it is
m x 1 vector, estimates of the variance function at each point of xgrid.
if type==2, it is
scalar (e.g., 1) and it has in this case no influence on the algorithm.

 msespan  optional scalar, number of points used to smooth the estimate of MSE.
If msespan is equal to 0 then the entire grid is used and all points
have the same weights. This yields an optimal global bandwidth. Setting
msespan to 0 is advocated for small sample size.

 nterms  optional scalar, number of terms in the bias model
The bias model is
E(fhat(xh)) = b_0+b_1 h^(p+1)+...+ b_nterms h^(p + nterms)

 bandspan  optional scalar, bandwidth for smoothing the bandwidth. By default, we do not
smooth the bandwidth.
If this smoothing is not desired then this parameter
should be negative (<= 0).

 varest  optional scalar, estimate of the variance function at each point of xgrid.

 Kernel  optional string. kernel function used in the local polynomial estimation.
The kernel functions available are Quartic ("Qua") (default),
Epanechnikov ("Epa") and Triangle ("Tri").

 J2  optional scalar, number of bandwidth points to the right used for
fitting a curve to estimate bias at one bandwidth.

 J1  optional scalar, number of bandwidth points to the left used for
fitting a curve to estimate bias at one bandwidth.

Output: 
 f  m x 1 vector representing the estimated OrderDerth derivative of E(yx) at each point
on xgrid. 
 band  m x 1 vector of bandwidths, entails the bandwidth
at each point of xgrid. 