bda (version 1.2.7-31)

tkde: Compute Transformation-Based Kernel Density Estimate

Description

Returns x and y coordinates of the kernel density estimate of the probability density of the data by transforming the data.

Usage

tkde(x, func.W, bw = "nrd0", adjust = 1, kernel = c("gaussian", 
    "epanechnikov", "rectangular", "triangular", "biweight", 
    "cosine", "optcosine"), window = kernel, width,
    n = 512, from, to, cut = 3, na.rm = FALSE, 
    ...)

Arguments

x

the data from which the estimate is to be computed.

func.W

The cumulative biasing functions.

bw

the smoothing bandwidth to be used. The kernels are scaled such that this is the standard deviation of the smoothing kernel. (Note this differs from the reference books cited below, and from S-PLUS.)

bw can also be a character s

adjust

the bandwidth used is actually adjust*bw. This makes it easy to specify values like half the default bandwidth.

kernel, window

a character string giving the smoothing kernel to be used. This must be one of "gaussian", "rectangular", "triangular", "epanechnikov", "biweight", "cosine" or

width

this exists for compatibility with S; if given, and bw is not, will set bw to width if this is a character string, or to a kernel-dependent multiple of width if this is numeric.

n

the number of equally spaced points at which the density is to be estimated. When n > 512, it is rounded up to a power of 2 during the calculations (as fft is used) and the final resul

from,to

the left and right-most points of the grid at which the density is to be estimated; the defaults are cut * bw outside of range(x).

cut

by default, the values of from and to are cut bandwidths beyond the extremes of the data. This allows the estimated density to drop to approximately zero at the extremes.

na.rm

logical; if TRUE, missing values are removed from x. If FALSE any missing values cause an error.

...

further arguments for (non-default) methods.

Value

• If give.Rkern is true, the number $R(K)$, otherwise an object with class "density" whose underlying structure is a list containing the following components.
• x
  the n coordinates of the points where the density is estimated.
• y
  the estimated density values. These will be non-negative, but can be zero.
• bw
  the bandwidth used.
• n
  the sample size after elimination of missing values.
• call
  the call which produced the result.
• data.name
  the deparsed name of the x argument.
• has.na
  logical, for compatibility (always FALSE).
• The print method reports summary values on the x and y components.

Details

The algorithm used in density.default disperses the mass of the empirical distribution function over a regular grid of at least 512 points and then uses the fast Fourier transform to convolve this approximation with a discretized version of the kernel and then uses linear approximation to evaluate the density at the specified points.

The statistical properties of a kernel are determined by $\sigma^2_K = \int t^2 K(t) dt$ which is always $= 1$ for our kernels (and hence the bandwidth bw is the standard deviation of the kernel) and $R(K) = \int K^2(t) dt$. MSE-equivalent bandwidths (for different kernels) are proportional to $\sigma_K R(K)$ which is scale invariant and for our kernels equal to $R(K)$. This value is returned when give.Rkern = TRUE. See the examples for using exact equivalent bandwidths.

Infinite values in x are assumed to correspond to a point mass at +/-Inf and the density estimate is of the sub-density on (-Inf, +Inf).

References

Becker, R. A., Chambers, J. M. and Wilks, A. R. (1988) The New S Language. Wadsworth & Brooks/Cole (for S version).

Scott, D. W. (1992) Multivariate Density Estimation. Theory, Practice and Visualization. New York: Wiley.

Sheather, S. J. and Jones M. C. (1991) A reliable data-based bandwidth selection method for kernel density estimation. J. Roy. Statist. Soc. B, 683--690.

Silverman, B. W. (1986) Density Estimation. London: Chapman and Hall.

Venables, W. N. and Ripley, B. D. (2002) Modern Applied Statistics with S. New York: Springer.

See Also

bw.nrd, plot.density, hist.

Examples

## biased sampling
 n = 1000; mu = 30.5; s = 8
 x = rnorm(n,mu,s); x = sort(x)
 fx = function(x,mu,s) dnorm(x,mu,s)*x/mu # length biasing
 x0 = seq(min(x)-s,max(x)+s,length=100)
 fn = fx(x0,mu,s)
 f0 = dnorm(x0, mu, s)

 plot(f0~x0, type='l',col=2, lwd=2)
 lines(fn~x0, lty=2, col=2, lwd=2)

 wx = bs.w.length(x); wx = wx/max(wx)
 sele = runif(length(x))<wx
 xw = x[sele]
 wx = 1/wx[sele]

 lines(density(xw), col=4, lwd=2)
 f1 = wkde(xw,w=wx)
 lines(f1, col=1, lwd=2)

 f2 = wkde(xw,w=wx, bandwidth='wmise')
 lines(f2, col=1, lty=2, lwd=2)

 f3 = wkde(xw,w=wx, bandwidth='awmise')
 lines(f3, col=1, lty=3, lwd=2)

 f4 = tkde(xw, func.W = bs.W.length)
 lines(f4, col=5, lty=1, lwd=2)

legend(max(x0), max(f0), xjust=1,yjust=1,
  legend=c("true", "biased","density(...)","f.rot","wmise","awmise","tkde"),
  lty=c(1,2,1,1,2,3,1), col=c(2,2,4,1,1,1,5),lwd=c(2,2,2,2,2,2,2))