Learn R Programming

mable (version 4.1.1)

mable.decon: Mable deconvolution with a known error density

Description

Maximum approximate Bernstein/Beta likelihood estimation in additive density deconvolution model with a known error density.

Usage

mable.decon(
  y,
  gn = NULL,
  ...,
  M,
  interval = c(0, 1),
  IC = c("none", "aic", "hqic", "all"),
  vanished = TRUE,
  controls = mable.ctrl(maxit.em = 1e+05, eps.em = 1e-05, maxit.nt = 100, eps.nt = 1e-10),
  progress = TRUE
)

Value

A mable class object with components, if \(g\) is known,

  • M the vector (m0, m1), where m1 is the last candidate degree when the search stoped

  • m the selected optimal degree m

  • p the estimate of p = (p_0, ..., p_m), the coefficients of Bernstein polynomial of degree m

  • lk log-likelihoods evaluated at \(m \in \{m_0, \ldots, m_1\}\)

  • lr likelihood ratios for change-points evaluated at \(m \in \{m_0+1, \ldots, m_1\}\)

  • convergence An integer code. 0 indicates an optimal degree is successfully selected in M. 1 indicates that the search stoped at m1.

  • ic a list containing the selected information criterion(s)

  • pval the p-values of the change-point tests for choosing optimal model degree

  • chpts the change-points chosen with the given candidate model degrees

if \(g\) is unknown,

  • M the 2 x 2 matrix with rows (m0, m1) and (k0,k1)

  • nu_aic the selected optimal degrees (m,k) using AIC method

  • p_aic the estimate of p = (p_0, ..., p_m), the coefficients of Bernstein polynomial model for \(f\) of degree m as in nu_aic

  • q_aic the estimate of q = (q_0, ..., q_k), the coefficients of Bernstein polynomial model for \(g\) of degree k as in nu_aic

  • nu_bic the selected optimal degrees (m,k) using BIC method

  • p_bic the estimate of p = (p_0, ..., p_m), the coefficients of Bernstein polynomial model for \(f\) of degree m as in nu_bic

  • q_bic the estimate of q = (q_0, ..., q_k), the coefficients of Bernstein polynomial model for \(g\) of degree k as in nu_bic

  • lk matrix of log-likelihoods evaluated at \(m \in \{m_0, \ldots, m_1\}\) and \(k \in \{k_0, \ldots, k_1\}\)

  • aic a matrix containing the Akaike information criterion(s) at \(m \in \{m_0, \ldots, m_1\}\) and \(k \in \{k_0, \ldots, k_1\}\)

  • bic a matrix containing the Bayesian information criterion(s) at \(m \in \{m_0, \ldots, m_1\}\) and \(k \in \{k_0, \ldots, k_1\}\)

Arguments

y

vector of observed data values

gn

error density function if known, default is NULL if unknown

...

additional arguments to be passed to gn

M

a vector (m0, m1) specifies the set of consective candidate model degrees, M = m0:m1. If gn is unknown then M a 2 x 2 matrix whose rows (m0,m1) and (k0,k1) specify lower and upper bounds for degrees m and k, respectively.

interval

a finite vector (a,b), the endpoints of supporting/truncation interval if gn is known. Otherwise, it is a 2 x 2 matrix whose rows (a,b) and (a1,b1) specify supporting/truncation intervals of X and \(\epsilon\), respectively. See Details.

IC

information criterion(s) in addition to Bayesian information criterion (BIC). Current choices are "aic" (Akaike information criterion) and/or "qhic" (Hannan–Quinn information criterion).

vanished

logical whether the unknown error density vanishes at both end-points of [a1,b1]

controls

Object of class mable.ctrl() specifying iteration limit and other control options. Default is mable.ctrl.

progress

if TRUE a text progressbar is displayed

Author

Zhong Guan <zguan@iu.edu>

Details

Consider the additive measurement error model \(Y = X + \epsilon\), where \(X\) has an unknown distribution \(F\) on a known support [a,b], \(\epsilon\) has a known or unknown distribution \(G\), and \(X\) and \(\epsilon\) are independent. We want to estimate density \(f = F'\) based on independent observations, \(y_i = x_i + \epsilon_i\), \(i = 1, \ldots, n\), of \(Y\). We approximate \(f\) by a Bernstein polynomial model on [a,b]. If \(g=G'\) is unknown on a known support [a1,b1], then we approximate \(g\) by a Bernstein polynomial model on [a1,b1], \(a1<0<b1\). We assume \(E(\epsilon)=0\). AIC and BIC methods are used to select model degrees (m,k).

References

Guan, Z., (2019) Fast Nonparametric Maximum Likelihood Density Deconvolution Using Bernstein Polynomials, Statistica Sinica, doi:10.5705/ss.202018.0173

Examples

Run this code
# \donttest{
 # A simulated normal dataset
 set.seed(123)
 mu<-1; sig<-2; a<-mu-sig*5; b<-mu+sig*5;
 gn<-function(x) dnorm(x, 0, 1)
 n<-50;
 x<-rnorm(n, mu, sig); e<-rnorm(n); y<-x+e;
 res<-mable.decon(y, gn, interval = c(a, b), M = c(5, 50))
 op<-par(mfrow = c(2, 2),lwd = 2)
 plot(res, which="likelihood")
 plot(res, which="change-point", lgd.x="topright")
 plot(xx<-seq(a, b, length=100), yy<-dnorm(xx, mu, sig), type="l", xlab="x",
     ylab="Density", ylim=c(0, max(yy)*1.1))
 plot(res, which="density", types=c(2,3), colors=c(2,3))
 # kernel density based on pure data
 lines(density(x), lty=4, col=4)
 legend("topright", bty="n", lty=1:4, col=1:4,
 c(expression(f), expression(hat(f)[cp]), expression(hat(f)[bic]), expression(tilde(f)[K])))
 plot(xx, yy<-pnorm(xx, mu, sig), type="l", xlab="x", ylab="Distribution Function")
 plot(res, which="cumulative",  types=c(2,3), colors=c(2,3))
 legend("bottomright", bty="n", lty=1:3, col=1:3,
     c(expression(F), expression(hat(F)[cp]), expression(hat(F)[bic])))
 par(op)
# }

Run the code above in your browser using DataLab