| iMqr {Mqrcm} | R Documentation |
M-Quantile Regression Coefficients Modeling
Description
This function implements Frumento and Salvati's (2020) method for M-quantile regression coefficients modeling (Mqrcm). M-quantile regression coefficients are described by parametric functions of the order of the quantile.
Usage
iMqr(formula, formula.p = ~ slp(p,3), weights, data, s,
psi = "Huber", plim = c(0,1), tol = 1e-6, maxit)
Arguments
formula |
a two-sided formula of the form |
formula.p |
a one-sided formula of the form |
weights |
an optional vector of weights to be used in the fitting process. The weights will always be normalized to sum to the sample size. This implies that, for example, using double weights will not halve the standard errors. |
data |
an optional data frame, list or environment containing the variables in |
s |
an optional 0/1 matrix that permits excluding some model coefficients (see ‘Examples’). |
psi |
a character string indicating the ‘psi’ function. Currently,
only ‘ |
plim |
the extremes of the estimation interval. You may want to model the M-quantile
regression coefficients in an interval, say, |
tol |
convergence criterion for numerical optimization. |
maxit |
maximum number of iterations. |
Details
A linear model is used to describe the p-th conditional M-quantile:
M(p | x) = \beta_0(p) + \beta_1(p)x_1 + \beta_2(p)x_2 + \ldots.
Assume that each M-quantile regression coefficient can be expressed as a parametric function of p of the form:
\beta(p | \theta) = \theta_{0} + \theta_1 b_1(p) + \theta_2 b_2(p) + \ldots
where b_1(p), b_2(p, \ldots) are known functions of p.
If q is the dimension of
x = (1, x_1, x_2, \ldots)
and k is that of
b(p) = (1, b_1(p), b_2(p), \ldots),
the entire M-conditional quantile function is described by a
q \times k matrix \theta of model parameters.
Users are required to specify two formulas: formula describes the regression model,
while formula.p identifies the 'basis' b(p).
By default, formula.p = ~ slp(p, k = 3), a 3rd-degree shifted
Legendre polynomial (see slp). Any user-defined function b(p, \ldots)
can be used, see ‘Examples’.
Estimation of \theta is carried out by minimizing an integrated loss function,
corresponding to the
integral, over p, of the loss function of standard M-quantile regression. This
motivates the acronym iMqr (integrated M-quantile regression). The scale parameter
sigma is estimated as the minimizer of the log-likelihood of a Generalized
Asymmetric Least Informative distribution (Bianchi et al 2017), and is “modeled”
as a piecewise-constant function of the order of the quantile.
Value
An object of class “iMqr”, a list containing the following items:
coefficients |
a matrix of estimated model parameters describing the fitted M-quantile function. |
plim |
a vector of two elements indicating the range of estimation. |
call |
the matched call. |
converged |
logical. The convergence status. |
n.it |
the number of iterations. |
obj.function |
the value of the minimized integrated loss function. |
s |
the used ‘s’ matrix. |
psi |
the used ‘ |
covar |
the estimated covariance matrix. |
mf |
the model frame used. |
PDF, CDF |
the fitted values of the conditional probability density function (PDF)
and cumulative distribution function (CDF). The CDF value should be interpreted as the order
of the M-quantile that corresponds to the observed |
Use summary.iMqr, plot.iMqr, and predict.iMqr
for summary information, plotting, and predictions from the fitted model.
The generic accessory functions coefficients, formula, terms,
model.matrix, vcov are available to extract information from the fitted model.
Author(s)
Paolo Frumento paolo.frumento@unipi.it
References
Frumento, P., Salvati, N. (2020). Parametric modeling of M-quantile regression coefficient functions with application to small area estimation, Journal of the Royal Statistical Society, Series A, 183(1), p. 229-250.
Bianchi, A., et al. (2018). Estimation and testing in M-quantile regression with application to small area estimation, International Statistical Review, 0(0), p. 1-30.
See Also
summary.iMqr, plot.iMqr, predict.iMqr,
for summary, plotting, and prediction, and plf and slp
that may be used to define b(p)
to be a piecewise linear function and a shifted Legendre polynomial basis, respectively.
Examples
##### Using simulated data in all examples
##### NOTE 1: the true quantile and M-quantile functions do not generally coincide
##### NOTE 2: the true M-quantile function is usually unknown, even with simulated data
##### Example 1
n <- 250
x <- runif(n)
y <- rnorm(n, 1 + x, 1 + x)
# true quantile function: Q(p | x) = beta0(p) + beta1(p)*x, with
# beta0(p) = beta1(p) = 1 + qnorm(p)
# fit the 'true' model: b(p) = (1 , qnorm(p))
m1 <- iMqr(y ~ x, formula.p = ~ I(qnorm(p)))
# the fitted M-quantile regression coefficient functions are
# beta0(p) = m1$coef[1,1] + m1$coef[1,2]*qnorm(p)
# beta1(p) = m1$coef[2,1] + m1$coef[2,2]*qnorm(p)
# a basis b(p) = (1, p), i.e., beta(p) is assumed to be a linear function of p
m2 <- iMqr(y ~ x, formula.p = ~ p)
# a 'rich' basis b(p) = (1, p, p^2, log(p), log(1 - p))
m3 <- iMqr(y ~ x, formula.p = ~ p + I(p^2) + I(log(p)) + I(log(1 - p)))
# 'slp' creates an orthogonal spline basis using shifted Legendre polynomials
m4 <- iMqr(y ~ x, formula.p = ~ slp(p, k = 3)) # note that this is the default
# 'plf' creates the basis of a piecewise linear function
m5 <- iMqr(y ~ x, formula.p = ~ plf(p, knots = c(0.1,0.9)))
summary(m1)
summary(m1, p = c(0.25,0.5,0.75))
par(mfrow = c(1,2)); plot(m1, ask = FALSE)
# see the documentation for 'summary.iMqr' and 'plot.iMqr'
##### Example 2 ### excluding coefficients
n <- 250
x <- runif(n)
qy <- function(p,x){(1 + qnorm(p)) + (1 + log(p))*x}
# true quantile function: Q(p | x) = beta0(p) + beta1(p)*x, with
# beta0(p) = 1 + qnorm(p)
# beta1(p) = 1 + log(p)
y <- qy(runif(n), x) # to generate y, plug uniform p in qy(p,x)
iMqr(y ~ x, formula.p = ~ I(qnorm(p)) + I(log(p)))
# I would like to exclude log(p) from beta0(p), and qnorm(p) from beta1(p)
# I set to 0 the corresponding entries of 's'
s <- rbind(c(1,1,0),c(1,0,1))
iMqr(y ~ x, formula.p = ~ I(qnorm(p)) + I(log(p)), s = s)