# Author: Aurore Delaigle
# Compute the measurement error version of the Nadaraya-Watson regression
# estimator
# Goal: estimate m where Y=m(X)+epsilon, and we observe data on (W, Y), where
# W=X+U.
# See Fan, J., and Truong, Y. K. (1993), Nonparametric Regression With
# Errors in Variables, The Annals of Statistics, 21, 1900-1925
# xx: vector of x-values where to compute the regression estimator
# W: vector of contaminated data W_1, ..., W_n
# Y: vector of data Y_1, ..., Y_n
# h: bandwidth
# rho: ridge parameter. See
# Delaigle, A. and Hall, P. (2008). Using SIMEX for smoothing-parameter choice
# in errors-in-variables problems. JASA, 103, 280-287
#
# -----------------------------------------------------------------------------
# WARNINGS:
# -----------------------------------------------------------------------------
#
# The phiK here must match the phiK used to compute the bandwidth (SIMEX or
# other).
#
# The estimator can also be computed using the Fast Fourier Transform, which
# is faster, but more complex.
# See Delaigle, A. and Gijbels, I. (2007). Frequent problems in calculating
# integrals and optimizing objective functions: a case study in density
# deconvolution. Statistics and Computing, 17, 349 - 355
# However if the grid of t-values is fine enough, the estimator can simply be
# computed like here without having problems with oscillations.
#
# -----------------------------------------------------------------------------
NWDecUknown <- function(xx, W, Y, phiU, h, rho, phiK, t, dt) {
# --------------------------------------------------------
# Preliminary calculations and initialisation of functions
# --------------------------------------------------------
W <- as.vector(W)
n <- length(W)
longt <- length(t)
# Compute the empirical characteristic function of W (times n) at t/h:
# \hat\phi_W(t/h)
OO <- t(outer(t/h, W))
csO <- cos(OO)
snO <- sin(OO)
rm(OO)
rehatphiW <- apply(csO, 2, sum)
imhatphiW <- apply(snO, 2, sum)
# Compute \sum_j Y_j e^{itW_j/h}
dim(Y) <- c(1, n)
renum <- Y %*% csO
imnum <- Y %*% snO
# Compute numerator and denominator of the estimator separately
# Numerator: real part of
# (2*pi*h)^(-1) \int e^{-itx/h} n^(-1)\sum_j Y_j e^{itW_j/h} \phi_K(t)/\phi_U(t/h) dt
#
# Denominator: real part of
# (2*pi*h)^(-1) \int e^{-itx/h} \hat\phi_W(t/h) \phi_K(t)/\phi_U(t/h) dt
xt <- outer(t / h, xx)
cxt <- cos(xt)
sxt <- sin(xt)
rm(xt)
phiUth <- phiU(t / h)
matphiKU <- phiK(t) / phiUth
dim(matphiKU) <- c(1, longt)
Den <- (rehatphiW * matphiKU) %*% cxt + (imhatphiW * matphiKU) %*% sxt
Num <- (renum * matphiKU) %*% cxt + (imnum * matphiKU) %*% sxt
Num <- Num * (dt / (n * h * 2 * pi))
Den <- Den * (dt / (n * h * 2 * pi))
# If denomintor is too small, replace it by the ridge rho
dd <- Den
dd[which(dd < rho)] <- rho
# Finally obtain the regression estimator
y <- Num / dd
as.vector(y)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.