modereg: Nonparametric Estimators for Nonparametric Mode Regression

Description

This function provides the nonparametric estimators (Zhou and Huang, 2016; Zhou and Huang, 2019) for nonparametric modal regression. The corresponding estimators in the absence of measurement error are also provided.

Usage

modereg(Y, W, bw, xgrid=NULL, sig=NULL, nstart=4, p.order=0, maxiter = 500, 
        tol=.Machine$double.eps^0.25, mesh=NULL, PLOT=FALSE, ...)

Arguments

an n by 1 response vector.

an n by 1 predictor vector.

bandwidth.

xgrid

the grid values to estimate the responses.

sig

standard deviation of the measurement error; sig=NULL returns the naive estimators ignoring measurement error.

nstart

the starting number of modes for each grid value.

p.order

the order of polynomial, up to 1; p.order=0 returns local constant estimators and p.order=1 returns local linear estimators.

maxiter

the maximum number of iterations performed for the mean shift algorithm if not converage.

tol

the deisered accurary (convergence tolerrance).

mesh

a matrix of initial mode points, where each row corresponds a mode in (x,y) coordinate; if mesh=NULL, it will be chosen automatically according to xgrid and nstart.

PLOT

a logical value indicating whether the estimated modes will be plotted.

...

further arguments to be passed to or from other methods.

Value

The results include the grid points xgrid for predictor, the number of modes for each grid x.num, the initial mesh points mesh, and corresponding fitted modes mode.

References

Zhou. H. and Huang, X. (2016). Nonparametric modal regression in the presence of measurement error. Electronic Journal of Statistics, 10: 3579-3620.

Zhou, H. and Huang, X. (2019). Bandwidth selection for nonparametric modal regression. Communications in Statistics - Simulation and Computation, 48(4): 968-984.

Examples

Run this code

# NOT RUN {
library(lpme)

rlaplace=function (use.n, location = 0, scale = 1) 
{
  location <- rep(location, length.out = use.n)
  scale <- rep(scale, length.out = use.n)
  rrrr <- runif(use.n)
  location - sign(rrrr - 0.5) * scale * 
    (log(2) + ifelse(rrrr < 0.5, log(rrrr), log1p(-rrrr)))
}

## sample size:
n =100;
## Function m(x) to estimate#
gofx1  = function(x){ (x+x^2) }
gofx2  = function(x){ (x+x^2)-6 }
xgrid	= seq(-2, 2, length.out=100); 
ngrid	= length(xgrid)

## Sample X
X = rnorm(n, 0, 1); sigma_x=1;
## Sample Y
Y = rep(0, n);
U = runif(n);
for(i in 1:n){
  if(U[i]<0.5){
    Y[i] = rnorm(1, gofx1(X[i]), 1);
  }else{
    Y[i] = rnorm(1, gofx2(X[i]), 1);
  }
}
## reliability ratio
lambda=0.9;
sigma_u  = sqrt(1/lambda-1)*sigma_x;
W=X+rlaplace(n,0,sigma_u/sqrt(2));

## mode estimates
hhxy = c(0.15, 1)
## Note you needs to use the following code to calculate bandwidth
## It is not run here due to the time constrain of runing examples. 
#hhxy = moderegbwSIMEX(Y, W, method="CV-density", p.order=0, 
#                       sig=sigma_u, B=5, length.h=10)$bw;
fit = modereg(Y, W, xgrid=xgrid, bw=hhxy, sig=sigma_u, p.order=0,
              PLOT=TRUE);

## Plot
plot(xgrid, gofx1(xgrid), "l", lwd="2", ylim=c(-9,7), xlim=c(-2,2));
lines(xgrid, gofx2(xgrid), "l", lwd="2");
points(rep(fit$xgrid,fit$x.num), fit$mode, col="3",lwd="2")
# }