drmdel (version 1.0)

densityDRM: Estimate the density of the populations under DRM

Description

Suppose we have m+1 samples, labeled as $0, \, 1, \, \ldots, \, m$, whose population distributions satisfy the density ratio model (DRM) (see drmdel for the definition of DRM). The densityDRM function estimates the density, $dF_k(x)$, of the k$\textsuperscript{th}$, $k=0, \, 1, \, \ldots, \, m$, population using a kernel density estimator with weights been the estimated $dF_k(x)$'s at the observed data points.

Usage

densityDRM(k, drmfit, interpolation=TRUE, adjust=FALSE,
           adj_factor=NULL, ...)

Arguments

k

the label of the population whose density is to be estimated, k = 0, 1, ..., m. It must be a single value, not a vector.

drmfit

a fitted DRM object (an output from the drmdel function). See drmdel for details.

interpolation

an argument passed to the function quantileDRM for calculating the default bandwidth for the kernel density estimator, if bandwidth (bw) is not specfied as a ... argument passing to the R fu

adjust

an argument passed to the function quantileDRM for calculating the default bandwidth for the kernel density estimator, if bandwidth (bw) is not specfied as a ... argument passing to R functi

adj_factor

an argument passed to the function quantileDRM for calculating the default bandwidth for the kernel density estimator, if bandwidth (bw) is not specfied as a ... argument passing to R functi

...

further arguments to be passed to the density function, which performs kernel density estimation in R. See help(density) for details. One should not pass "x" and "weights" arguments to the density functio

Value

• An output from the density function, usually an object with class '"density"'. See help(density) for details.

Details

Note that the default bandwidth for this density estimator is set as that described in Chen and Liu (2013): $$1.06 n_{total}^{-1/5} \min(\sigma_k, \, \mbox{IQR}_k/1.34)$$, where $n_{total}$ is the total sample size, and $sigma_k$ and $IQR_k$ are the standard deviation and inter-quartile range of the estimated CDF $F_k$, respectively.

If bandwidth (bw) is not specfied as a ... argument passing to the R functin density, the default bandwidth (as described above) will be calculated. That calculation involves the estimation of population quartiles. In this situation, the arguments "interpolation", "adjust", and "adj_factor" are passed to quantileDRM for quartile estimations. See quantileDRM for details.

References

K. Fokianos (2004), Merging information for semiparametric density estimation. Journal of the Royal Statistical Society, Series B (Statistical Methodology), 66(4):941-958.

J. Chen and Y. Liu (2013), Quantile and quantile-function estimations under density ratio model. To appear in The Annals of Statistics, 2013.

Examples

# Data generation
set.seed(25)
n_samples <- c(100, 200, 180, 150, 175)  # sample sizes
x0 <- rgamma(n_samples[1], shape=5, rate=1.8)
x1 <- rgamma(n_samples[2], shape=12, rate=1.2)
x2 <- rgamma(n_samples[3], shape=12, rate=1.2)
x3 <- rgamma(n_samples[4], shape=18, rate=5)
x4 <- rgamma(n_samples[5], shape=25, rate=2.6)
x <- c(x0, x1, x2, x3, x4)

# Fit a DRM with the basis function q(x) = (x, log(abs(x))), which
# is the basis function for gamma family. This basis function is
# the built-in basis function 6.
drmfit <- drmdel(x=x, n_samples=n_samples, basis_func=6)

# Estimate the density of population 3 under the DRM
dens_pop3 <- densityDRM(k=3, drmfit=drmfit)

# Plot the estimated density
plot(dens_pop3, main=bquote(F[3]), ylim=range(c(0, 0.5)))

# Add the empirical kernel density estimation curve of F_3 based
# on the third sample on the above density plot
lines(density(x3), col="blue", lty="28F8")

# Add the true density curve of F_3 on the above density plot
lines(seq(min(dens_pop3$y), max(dens_pop3$x), 0.01),
      dgamma(seq(min(dens_pop3$y), max(dens_pop3$x), 0.01), 18,
             5),
      type="l", col="red", lty="dotted")

legend(9, 0.5,
       legend=c("DRM density estimator",
                "Empirical kernel density estimator",
                "True density"),
       col=c("black", "blue", "red"),
       lty=c("solid", "28F8", "dotted"))