qde: Quantiles based Density Estimation

In qrjoint: Joint Estimation in Linear Quantile Regression

Quantiles based Density Estimation

Description

Provides a semiparametric estimation of the quantiles for independented univariate data with possible right censoring. This is same as estimating the intercept function within a joint linear quantile regression model with no predictors.

Usage

qde(y, nsamp = 1e3, thin = 10, cens = NULL,
    wt = NULL, incr = 0.01, par = "prior", nknots = 6,
    hyper = list(sig = c(.1,.1), lam = c(6,4), kap = c(0.1,0.1,1)),
    prox.range = c(.2,.95), acpt.target = 0.15, ref.size = 3,
    blocking = "single", temp = 1, expo = 2, blocks.mu, blocks.S,
    fix.nu = FALSE, fbase = c("t", "logistic", "unif"), verbose = TRUE)

## S3 method for class 'qde'
update(object, nadd, append = TRUE, ...)

Arguments

y	numeric vector of response data.
nsamp	number of posterior samples to be saved; defaults to 1000.
thin	thinning rate for the Markov chain sampler – one posterior sample is saved per thin iterations. Defaults to 10. The Markov chain sampler runs for a total of nsamp * thin many iterations.
cens	censoring status of response. Must be a vector of length length(y), with 0 indicating no censoring, 1 indicating right censoring, and 2 indicating left censoring. If not supplied, defaults to all zeros.
wt	weights attached to the observation units, expected to be non-negative numbers, and defaults to a vector of ones if not otherwise supplied.
incr	tau grid increment. Defaults to 0.01.
par	character string indicating how the sampler is to be initialized. Only two options are currently supported: "prior" to initialize at a random draw from the prior; "RQ" to initialize at a model space approximation of the estimates from rq.
nknots	number of knots to be used for low rank approximation of the Gaussian process priors. Defaults to 6.
hyper	hyperparameters of the prior distribution. Must be a list with some of all of the following fields: sig: a two vector giving the parameters of the inverse-gamma distribution on sigma-square that is used when shrink=TRUE, lam: a two vector giving the parameters of the beta distribution on proximity = \exp(-0.01* \lambda^2), and kap: a vector to be coerced into a 3 * nkap matrix, with nkap being the number of components in the mixture of gamma prior on kappa, and each column of the matrix gives the shape, rate and mixing weight of a component.
prox.range	for specifying the range of length-scale parameter of the Gaussian process prior.
acpt.target	target acceptance rate of the adaptive Metropolis sampler; defaults to 0.15
ref.size	adaptation rate of the adaptive Metropolis sampler. The proposal density is updated once every ref.size iterations. Could be a single number or a vector of length same as the number of blocks.
blocking	type of blocking to be applied. Either a character string specifying one to be chosen from the supplied menu (see Details), or a list giving user specified blocks. In the latter case, each element of the list is a logical vector of length equal to the total number of model parameters, which equals nknots+3 indicating which model parameters belong to the block.
temp	temperature of the log-likelihood function. The log-likelihood function is raised to the power of temp. Defaults to 1.
expo	the exponent to be used in the covariance kernel of the Gaussian process priors. Defaults to 2, giving the standard squared-exponential covariance kernel.
blocks.mu	initial block specific means in the form of a list. If left unspecified then will be automatically generated as a list of vectors of zeros of appropriate lengths matching the corresponding block sizes.
blocks.S	initial block specific covariance matrices in the form of a list. If left unspecified then will be automatically generated as a list of identity matrices of appropriate dimensions matching the corresponding block sizes. When blocking is chosen as one of the menu items of the form "std*", known prior covariance information and estimated variance matrices from rq are used.
fix.nu	either the logical FALSE indicating that nu should be learned, or a positive real number giving the fixed value of nu, which is then excluded from MCMC updates
fbase	either "t" (default), "logistic" or "unif" to indicate what base distribution is to be used.
verbose	logical indicating whether MCMC progress should be printed, defaults to TRUE
object	a fitted model of the class 'qde'.
nadd	number of additional MCMC samples.
append	logical indicating whether new samples should be appended to old ones. If FALSE then old samples are discarded.
...	no additional arguments are allowed

Details

The model assumes the quantile function of the data is given by: Q(t) = gamma_0 + sigma * (Qb(zeta(t)|nu) - Qb(zeta(0,5)|nu)), Qb(.|nu) is a parametric quantile function with unknown parameter nu, gamma_0 is the unknown median, sigma is an unknown saling factor, and, zeta is an unknown distortion of the unit interval. The distortion zeta is modeled nonparametrically through a logistic Gaussian process prior, other parameters are given diffuse priors.

In running the MCMC, the following menu choices are available for blocking the parameter vector. For this special case p = ncol(X) = 0, some of the menu choices are actually the same, in particular, "std0" is same as "single", "std1" is same as "single2", and, "std2" is same as "single3".

"single": a single block containing all parameters

"single2": one block containing all parameters and an additional block containing only (gamma[0], gamma, sigma, nu)

"single3": like "single2", but the second block is split into two further blocks, one with (\gamma_0, \gamma), the other with (\sigma, \nu)

"std0": Same as "single".

"std1": Same as "single2".

"std2": Same as "single3".

"std3": total 3 blocks. First block is W_{*0}, last two are (\gamma_0, \gamma) and (\sigma, \nu)

"std4": total 3 blocks. First block is (W_{*0}, \gamma_0), last two are (\gamma_0, \gamma) and (\sigma, \nu)

"std5": total 4 blocks. First three are same as "std4" and one additional block containing all parameters.

Value

qde(y, ...) returns a ‘qde’ class object to be used by coef and summary.

Returned object is a list containing the following variables.

par	latest draw of the parameter vector
y	response vector
cens	censoring status vector
wt	vector of observation weights
hyper	completed list of hyper-parameters
dim	model dimension vector of the form c(n, p, length of tau grid, position of \tau_0 on the grid, nknots, length of lambda grid, nkap, total number of MCMC iterations, thin, nsamp)
gridmats	details of covariance matrix factors etc, intended for internal use.
tau.g	the tau grid
muV	list of means for parameter blocks
SV	list of covariance matrices for parameter blocks
blocks	list of blocks
blocks.size	vector of block lengths
dmcmcpar	numeric vector containing details of adaptive MCMC runs, equals c(temp, decay rate of adaptation, vector of target acceptance rates for the blocks, vector of increment scales used in adaptation). Intended strictly for internal use.
imcmcpar	numeric vector containing details of adaptive MCMC runs, equals c(number of parameter blocks, ref.size, indicator on whether details are to be printed during MCMC progress, rate of details printing, a vector of counters needed for printing). Intended strictly for internal use.
parsamp	a long vector containing the parameter draws. Could be coerced into a matrix of dim npar * nsamp. Intended primarily for use by summary and coef.
acptsamp	a long vector containing rates of acceptance statistics for parameter blocks. Could be coerced into a matrix of dim nblocks * nsamp. Not very informative, because thinning times and adaptation times may not be exactly synced.
lpsamp	vector of log posterior values for the saved MCMC draws.
fbase.choice	integer 1 for "t", 2 for "logistic" and 3 for "unif" base.
prox	vector of proximity (exp(-0.01*lambda^2)) grid values
reg.ix	positions of the regular tau grid on the expanded tail-appended grid
runtime	run time of the MCMC

References

Yang, Y. and Tokdar, S.T., 2017. Joint estimation of quantile planes over arbitrary predictor spaces. Journal of the American Statistical Association, 112(519), pp.1107-1120.

See Also

summary.qde, coef.qde and predict.qde. Also see qrjoint for regression model fitting in presence of covariates.

Examples

## Plasma data analysis

data(plasma)
Y <- plasma$BetaPlasma

# model fitting with 100 posterior samples from 200 iterations (thin = 2)
# this is of course for illustration, for practical model fitting you
# ought to try at least something like nsamp = 500, thin = 20
fit.qde <- qde(Y, nsamp = 100, thin = 2)
summary(fit.qde, more = TRUE)
pred <- predict(fit.qde)
hist(Y, freq = FALSE, col = "gray", border = "white", ylim = c(0, max(pred$fest)))
lines(pred$y, pred$fest[,2])
lines(pred$y, pred$fest[,1], lty = 2)
lines(pred$y, pred$fest[,3], lty = 2)

qrjoint documentation built on April 6, 2023, 1:07 a.m.