| Title: | Valid Inference on Multiple Quantile Regressions |
|---|---|
| Description: | The approach is based on the closed testing procedure to control familywise error rate in a strong sense. The local tests implemented are Wald-type and rank-score. The method is described in De Santis, et al., (2026), <doi:10.48550/arXiv.2511.07999>. |
| Authors: | Angela Andreella [aut, cre], Anna Vesely [aut], Riccardo De Santis [aut] |
| Maintainer: | Angela Andreella <[email protected]> |
| License: | GPL (>= 3) |
| Version: | 0.2.1 |
| Built: | 2026-05-10 07:57:02 UTC |
| Source: | https://github.com/angeella/quasar |
The Imhof method for x=0 and central variables only.
.pImhof(lams, eps = c(1e-04, 1e-04), x).pImhof(lams, eps = c(1e-04, 1e-04), x)
lams |
eigenvalues |
eps |
tolerance vector (one for integrate, one for quadgk) |
x |
observed stat test |
Applies the closed testing procedure to strongly control the familywise error rate (FWER) when testing the effect of a covariate of interest across multiple quantile regression models.
closedTesting(mod, X, tau = NULL, test = "rank-score", ...)closedTesting(mod, X, tau = NULL, test = "rank-score", ...)
mod |
An object of class |
X |
A string indicating the covariate of interest. |
tau |
A numeric vector of quantiles of interest used in |
test |
Character. Type of test to be used. Options are
|
... |
This procedure requires that the covariate of interest X is either numeric
or, if categorical, has at most two levels. Multilevel categorical covariates
are not supported and will trigger an error.
The weighting matrix used in the multivariate rank test is the identity matrix, i.e., it is currently the only one implemented explicitly as a shortcut within the closed testing procedure.
An object of class quasar containing:
Quantile: quantile level
Coefficient: estimated coefficient
Statistic: test statistic
p.value: raw -value
p.value.adjusted: adjusted -value from the closed testing procedure
Angela Andreella
Marcus, R., Eric, P., & Gabriel, K. R. (1976). On closed testing procedures with special reference to ordered analysis of variance. Biometrika, 63(3), 655–660.
Goeman, J. J., Hemerik, J., & Solari, A. (2021). Only closed testing procedures are admissible for controlling false discovery proportions. The Annals of Statistics, 49(2), 1218–1238.
# Simulate data set.seed(1234) D <- simulateData(n = 100, gamma = 0.5, sigma.y = "1 + 2 * pmax(X, 0)") # Quantile regressions at different levels tau <- c(0.1, 0.25, 0.5, 0.75, 0.9) mod <- quantreg::rq(y ~ X + Z1, tau = tau, data=D) # Closed testing res <- closedTesting(mod, X = "X") res # Summary and plot summary(res, alpha = 0.1) plot(res, alpha = 0.1, legend.position = "bottomright")# Simulate data set.seed(1234) D <- simulateData(n = 100, gamma = 0.5, sigma.y = "1 + 2 * pmax(X, 0)") # Quantile regressions at different levels tau <- c(0.1, 0.25, 0.5, 0.75, 0.9) mod <- quantreg::rq(y ~ X + Z1, tau = tau, data=D) # Closed testing res <- closedTesting(mod, X = "X") res # Summary and plot summary(res, alpha = 0.1) plot(res, alpha = 0.1, legend.position = "bottomright")
Produces a plot of a quasar object,
typically returned by the closedTesting function.
It shows the estimated coefficients by quantile level, highlighting
statistically significant coefficients based on adjusted p-values.
## S3 method for class 'quasar' plot( x, alpha = 0.05, legend.position = "topright", main = NULL, xlab = "Quantile level", ylab = "Coefficient", col.line = "darkgrey", col.sig = "darkred", col.nonsig = "darkgrey", pch.sig = 19, pch.nonsig = 17, show.legend = TRUE, ... )## S3 method for class 'quasar' plot( x, alpha = 0.05, legend.position = "topright", main = NULL, xlab = "Quantile level", ylab = "Coefficient", col.line = "darkgrey", col.sig = "darkred", col.nonsig = "darkgrey", pch.sig = 19, pch.nonsig = 17, show.legend = TRUE, ... )
x |
An object of class |
alpha |
Significance level. |
legend.position |
Position of the legend. |
main |
Main plot title. |
xlab |
Label for the x-axis. |
ylab |
Label for the y-axis. |
col.line |
Color of the connecting line. |
col.sig |
Color for significant points. |
col.nonsig |
Color for non-significant points. |
pch.sig |
Point character for significant points. |
pch.nonsig |
Point character for non-significant points. |
show.legend |
Logical; whether to display a legend. |
... |
Additional graphical parameters passed to |
A base R plot.
Anna Vesely
These methods provide basic information about objects of class
quasar, typically returned by the closedTesting function.
## S3 method for class 'quasar' print(x, ...) ## S3 method for class 'quasar' summary(object, ..., alpha = 0.05)## S3 method for class 'quasar' print(x, ...) ## S3 method for class 'quasar' summary(object, ..., alpha = 0.05)
x, object
|
An object of class |
... |
Additional arguments passed to other methods. |
alpha |
Significance level. |
The input object invisibly.
Anna Vesely
Performs the rank-score test for the covariate of interest
X, at the quantiles defined in tau, using a fitted quantile
regression model. The test evaluates the null hypothesis that the coefficient of X is equal to zero
against a two-sided alternative, at each specified quantile level.
Testing equality to a non-zero value is not yet implemented.
rankTest(mod, X, tau = NULL, full = FALSE, h = NULL, alpha = 0.05, eps = c(1e-04,1e-04), B = "identity", error.distr = NULL, error.par = NULL)rankTest(mod, X, tau = NULL, full = FALSE, h = NULL, alpha = 0.05, eps = c(1e-04,1e-04), B = "identity", error.distr = NULL, error.par = NULL)
mod |
An object of class |
X |
A string indicating the coefficient of interest, i.e. the name of a column in the model matrix |
tau |
A numeric vector of quantiles of interest used in |
full |
Logical. If |
h |
Character string specifying the bandwidth selection method for
sparsity estimation. One of |
alpha |
A numeric value used for bandwidth estimation. Following Koenker (2005), it is typically set equal to the desired significance level. |
eps |
Numeric vector of length 2 specifying the relative accuracies
requested for approximating the distribution of the test statistic using
Imhof's (1961) method, "Computing the Distribution of Quadratic Forms in
Normal Variables". The first component controls the accuracy of the initial
(more accurate) numerical integration (using |
B |
Weight specification used in the computation of the test statistic.
One of |
error.distr |
A character string specifying the assumed distribution of the
regression errors, used only when |
error.par |
A named list of parameters associated with
|
This procedure requires that the covariate of interest X is either numeric
or, if categorical, has at most two levels. Multilevel categorical covariates
are not supported and will trigger an error.
If full = TRUE, B is ignored and the identity weight is used.
If B is user-defined (a matrix), it must be square, symmetric, positive
definite, and numerically invertible.
A data.frame containing:
Quantiles.Set: quantile levels
Statistic: rank-score test statistic
p.value: corresponding test-specific (unadjusted) -value
Angela Andreella
Koenker, R. (2005). Quantile Regression. Cambridge University Press.
Imhof, J. P. (1961). Computing the distribution of quadratic forms in normal variables. Biometrika, 48(3–4), 419–426.
De Santis, R., Veselý, A., and Andreella, A. (2026). Inference on multiple quantiles in regression models by a rank-score approach. arXiv preprint <doi:10.48550/arXiv.2511.07999>.
set.seed(1234) D <- simulateData(n = 100, gamma = 0.5, sigma.y = "1 + 2 * pmax(X, 0)") #Quantile regressions at different levels tau <- c(0.1, 0.25, 0.5, 0.75, 0.9) mod <- quantreg::rq(y ~ X + Z1, tau = tau, data=D) # Rank test rankTest(mod, X = "X")set.seed(1234) D <- simulateData(n = 100, gamma = 0.5, sigma.y = "1 + 2 * pmax(X, 0)") #Quantile regressions at different levels tau <- c(0.1, 0.25, 0.5, 0.75, 0.9) mod <- quantreg::rq(y ~ X + Z1, tau = tau, data=D) # Rank test rankTest(mod, X = "X")
Simulates a main covariate X,
a vector of additional covariates Z, and a response y drawn from
the chosen distribution.
simulateData(n, beta = 0, gamma = 0, mu = 0, Sigma = NULL, sigma.y = 1, distribution = "normal", df = 5, xi = -1.453, omega = 2, alpha = 2.2, seed = NULL)simulateData(n, beta = 0, gamma = 0, mu = 0, Sigma = NULL, sigma.y = 1, distribution = "normal", df = 5, xi = -1.453, omega = 2, alpha = 2.2, seed = NULL)
n |
Integer. Number of observations. |
beta |
Numeric scalar. Effect of |
gamma |
Numeric vector. Effects of |
mu |
Numeric scalar. Intercept. |
Sigma |
Numeric |
sigma.y |
Either a numeric scalar or a one-sided expression/string (e.g., |
distribution |
Character. One of |
df |
Numeric scalar > 0. Degrees of freedom for t-distribution. |
xi |
Numeric scalar. Location parameter for the skew-normal distribution. In particular, this will be
|
omega |
Numeric scalar > 0. Scale parameter for the skew-normal distribution. In particular, this will be
|
alpha |
Numeric scalar. Slant parameter for the skew-normal distribution. Default 2.2. |
seed |
Numeric scalar > 0. Seed for random number generator. |
The response is generated as y = mu + X * beta + Z %*% gamma + error.
The error term can be drawn from a normal distribution, scaled Student-t with df degrees of freedom,
or a skew-normal. Its standard deviation is defined by sigma.y:
if numeric, a fixed scale is used; if a character expression,
the scale can vary with X and/or Z1.
A data.frame with columns y, X, and Z1, ..., Zk.
Angela Andreella
set.seed(1) p <- 3 Sigma <- diag(p) # Normal dat_n <- simulateData(n = 200, beta = 0.5, gamma = c(0.2,-0.1), sigma.y = 0.5, distribution = "normal") # Student-t dat_t0 <- simulateData(n = 200, beta = 0.5, gamma = c(0.2,-0.1), sigma.y = 0.5, distribution = "t", df = 7) # Skew-normal dat_sn <- simulateData(n = 200, beta = 0.5, gamma = c(0.2,-0.1), sigma.y = "abs(Z1) + 1", distribution = "skew-normal")set.seed(1) p <- 3 Sigma <- diag(p) # Normal dat_n <- simulateData(n = 200, beta = 0.5, gamma = c(0.2,-0.1), sigma.y = 0.5, distribution = "normal") # Student-t dat_t0 <- simulateData(n = 200, beta = 0.5, gamma = c(0.2,-0.1), sigma.y = 0.5, distribution = "t", df = 7) # Skew-normal dat_sn <- simulateData(n = 200, beta = 0.5, gamma = c(0.2,-0.1), sigma.y = "abs(Z1) + 1", distribution = "skew-normal")
Performs the Wald-type test for the covariate of interest
X, at the quantiles defined in tau, using a fitted quantile
regression model. The test evaluates the null hypothesis that the coefficient of X is equal to
a given value beta against a two-sided alternative, at each specified quantile level.
waldTest(mod, X, tau = NULL, full = FALSE, h = NULL, beta = 0, alpha = 0.05)waldTest(mod, X, tau = NULL, full = FALSE, h = NULL, beta = 0, alpha = 0.05)
mod |
An object of class |
X |
A string indicating the covariate of interest. |
tau |
A numeric vector of quantiles of interest used in |
full |
Logical. If |
h |
A numeric value for the bandwidth. |
beta |
Numeric value of the parameter of interest under the null hypothesis. |
alpha |
A numeric value used for bandwidth estimation. Following Koenker (2005), it is typically set equal to the desired significance level. |
This procedure requires that the covariate of interest X is either numeric
or, if categorical, has at most two levels. Multilevel categorical covariates
are not supported and will trigger an error.
A data.frame containing:
Quantiles.Set: quantile levels
Statistic: Wald-type test statistic
p.value: corresponding test-specific (unadjusted) -value
Angela Andreella
Koenker, R. (2005). Quantile Regression. Cambridge University Press.
set.seed(1234) D <- simulateData(n = 100, gamma = 0.5, sigma.y = "1 + 2 * pmax(X, 0)") #Quantile regressions at different levels tau <- c(0.1, 0.25, 0.5, 0.75, 0.9) mod <- quantreg::rq(y ~ X + Z1, tau = tau, data=D) # Wald test waldTest(mod, X = "X")set.seed(1234) D <- simulateData(n = 100, gamma = 0.5, sigma.y = "1 + 2 * pmax(X, 0)") #Quantile regressions at different levels tau <- c(0.1, 0.25, 0.5, 0.75, 0.9) mod <- quantreg::rq(y ~ X + Z1, tau = tau, data=D) # Wald test waldTest(mod, X = "X")