# R/meanreg.R In lpme: Nonparametric Estimation of Measurement Error Models

```"meanreg" <- function(Y, W, bw, xgrid=NULL, method="HZ", sig=NULL, error="laplace", FT_fu){
#########################################################################################
# data structure
#########################################################################################
if(is.null(xgrid)){
xgrid = seq(quantile(W, probs=0.025), quantile(W,probs=0.975), length.out = 100);
}
xmin  = min(xgrid);
xmax  = max(xgrid);
m  = 2^16; m_mid = m/2+1;
beta  = sqrt((2*pi)/m);
beta2  = (2*pi)/(m*beta);
input  = seq(-m*beta/2, m*beta/2-beta, by=beta)
output= seq(-pi/beta, pi/beta-beta2, by=beta2)
nlower= round(xmin/beta+m/2+1); nupper = round(xmax/beta+m/2+1); xindex = nlower:nupper;
mconst= (-1)^(0:(m-1)); xindexC = xindex-1;
x  = input[xindex];
wide	= beta;
nx	= length(x);
## Kernel for which CF is (1-t^2)^8 with Normal errors
FKsup  = function(t) { ifelse( (t<=1 & t>=-1), (1-t^2)^8, 0) }
FKoutput = rep(0, m);
FKoutput[m_mid] = FKsup(output[m_mid]);
i=1; indicator1 =1; indicator2 =1;
while ( ( abs(indicator1)>1e-30 | abs(indicator2)>1e-30 ) & (i<(m/2)) ) {
indicator1 = FKsup(output[m_mid-i]);
indicator2 = FKsup(output[m_mid+i]);
FKoutput[m_mid-i] = indicator1;
FKoutput[m_mid+i] = indicator2;
i=i+1;
}
## inverse FFT to get K
Kinput  = Re( (-1)^(0:(m-1))/beta*fft( (-1)^(0:(m-1))*FKoutput)/m );
dt = 0.0001;
tt = seq(-1,1,dt);

#########################################################################################
# calling the c++ code and # output
#########################################################################################
if(method=="HZ"){
if(error=="laplace"){
if(is.null(sig)) stop("please specify sig if non naive method is used")
foo <- .Call("fitnewLap",
x_ = x,
input_ = input,
output_ = output,
beta_ = beta,
beta2_ = beta2,
mconst_ = mconst,
Kinput_ = Kinput,
W_ = W,
Y_ = Y,
sigU_ = sig,
h_ = bw,
PACKAGE = "lpme");
}else if(error=="normal") {
if(is.null(sig)) stop("please specify sig if non naive method is used")
foo <- .Call("fitnewGau",
x_ = x,
input_ = input,
output_ = output,
beta_ = beta,
beta2_ = beta2,
mconst_ = mconst,
Kinput_ = Kinput,
W_ = W,
Y_ = Y,
sigU_ = sig,
h_ = bw,
PACKAGE = "lpme");
}else{
if(!is.function(FT_fu)) stop("pleaes specify the characteristic function of U as a function");
FfU = FT_fu(output);
foo <- .Call("fitNEW",
x_ = x,
input_ = input,
output_ = output,
beta_ = beta,
beta2_ = beta2,
mconst_ = mconst,
Kinput_ = Kinput,
W_ = W,
Y_ = Y,
FfU_ = FfU,
h_ = bw,
PACKAGE = "lpme");
}
} else if (method=="DFC") {
if(is.null(sig)) stop("please specify sig if non naive method is used")
if(error=="laplace"){
foo <- .Call("fitjasaLap",
x_ = x,
h_ = bw,
W_ = W,
Y_ = Y,
sigU_ = sig,
dt_ = dt,
t_ = tt,
PACKAGE = "lpme");
} else {
foo <- .Call("fitjasaGau",
x_ = x,
h_ = bw,
W_ = W,
Y_ = Y,
sigU_ = sig,
dt_ = dt,
t_ = tt,
PACKAGE = "lpme");
}
} else {
foo <- .Call("fitlocpoly",
x_ = x,
beta_ = beta,
Kinput_ = Kinput,
W_ = W,
Y_ = Y,
h_ = bw,
PACKAGE = "lpme");
}

#########################################################################################
# save state
#########################################################################################
model.name <- "local polynomial estimation";
ngrid = length(xgrid);
yhat = rep(0, ngrid);
for(i in 1:ngrid){
yhat[i] = (foo\$ghat)[which.min(abs(x-xgrid[i]))]
}

#### Save to a list
output <- list(modelname=model.name,
xgrid = xgrid,
yhat = yhat);
class(output) <- c("meanreg")
output
}
```

## Try the lpme package in your browser

Any scripts or data that you put into this service are public.

lpme documentation built on July 14, 2021, 1:06 a.m.