*! version 1.0.1 14apr2020

Title

qrprocess -- Quantile regression: fast algorithms, pointwise and uniform inference

Syntax

qrprocess depvar [indepvars] [if] [in]
[weight]
[, options]

options

Quantiles to be estimated
quantile(numlist) specifies the quantile(s) to be estimated; default value is 0.5 (the median) except if the option functional is activated.

Algorithmic method
method([mtype], [mopts]) specifies the algorithmic method used to compute the quantile regression estimates as well as the parameters governing the algorithm.

Inference
vce([vtype], [vopts]) specifies the technique used to estimate standard errors (variance-covariance matrix of the estimates).

functional activates the tools for functional inference (uniform confidence bands, test of functional hypotheses).

Reporting
level(#) sets confidence level; default is level(95).

noprint suppresses display of the results.

mtype

qreg selects the built-in Stata algorithm used by the official command qreg; the default for sample sizes below 10,000 observations when less than 10 quantile regressions must be estimated.

fn selects a basic implementation of the Frisch-Newton interior point method suggested by Portnoy and Koenker (1997).

pqreg selects the built-in Stata algorithm after the preprocessing step suggested by Portnoy and Koenker (1997).

pfn selects the Frisch--Newton interior point method after the preprocessing step suggested by Portnoy and Koenker (1997); the default for sample sizes above 10,000 observations when less than 10 quantile regressions must be estimated.

proqreg selects the built-in Stata algorithm after the preprocessing step based on the previously estimated quantile regression (preprocessing for the quantile regression process).

profn selects the Frisch--Newton interior point method after the preprocessing step based on the previously estimated quantile regression (preprocessing for the quantile regression process); the default for sample sizes below 50'000 when at least 10 quantile regressions are estimated.

1step selects the one-step estimator that uses a first-order approximation from one quantile regression to the next; the default for sample sizes above 50'000 when at least 10 quantile regressions are estimated.

mopts

beta(#) technical step length parameter used by the fn method; default value is 0.9995.

small(#) tolerance parameter for convergence used by the fn method; default value is 0.00001.

max_it(#) maximum number of iterations for fn method; default value is 100.

m(#) adjustement parameter used by the preprocessing algorithms; default value is 3.

error_tol(#) integer determining the number of allowed mispredicted signs of residuals for the methods involving preprocessing and for the bootstrap. The default value is 0.

first(mtype) selects the algorithm for the first quantile regression when the method is proqreg, profn or 1step.

vtype

kernel selects the kernel estimates of the VCE suggested by Powell (1991); the default when functional tests and confidence bands are not requested.

nid selects the estimate of the VCE suggested by Hendricks and Koenker (1992).

bootstrap selects a bootstrap method; the default when functional tests and confidence bands are requested.

multiplier multiplier bootstrap estimates of the VCE.

iid selects an estimate of the VCE that assumes iid standards errors. Not recommended!

novar does not compute the VCE.

vopts

bofinger selects Bofinger's bandwidth. By default, Hall-Sheather's bandwidth is used.

misspecification selects the Angrist, Chernozhukov and Fernández-Val (2006) estimator of the variance, which is robust against misspecification.

cluster(varname) cluster identification variable.

strata(varname) strata identification variable.

bmethod(bmethod) bootstrap method.

reps(#) number of bootstrap replications; default value is 100.

noreplacement subsampling only: performs subsampling without replacement; default is subsampling with replacement.

subsize(#) subsampling only: size of the subsamples.

functional requests functional inference results.

bmethod

empirical uses the empirical bootstrap; the default with vce(bootstrap).

subsampling uses subsampling.

weighted uses the weighted bootstrap with standard exponential weights; can only be used with vce(bootstrap).

wild uses the wild bootstrap; the default with vce(multiplier); can only be used with vce(multiplier).

Gaussian uses Gaussian weights; can only be used with vce(multiplier).

weighted uses centered standard exponential weights; can only be used with vce(multiplier).

by is allowed with qrprocess.

qrprocess allows fweights, iweights, and pweights; see weight.

See qrprocess postestimation for features available after estimation.

Description

qrprocess fits linear quantile regression models. This commands implements several algorithms, including new algorithms that are much quicker than
the built-in Stata commands, especially for big data or when a large number of regressions or bootstrap replications must be estimated. qrprocess also provides
analytical estimates of the variance-covariance matrix of the coefficients for several regressions allowing for weights, clustering and stratification. In addition
to traditional pointwise confidence intervals, these commands also provide functional confidence bands and tests of functional hypotheses.
depvar is the dependent variable. indepvar the regressor(s). The quantile(s) to be estimated are determined by the option quantile.
The algorithmic method used to obtain the estimates is determined by the option method.
The pointwise standard errors and confidence intervals are calulated according to the option vce.
Functional tests and confidence bands are provided if the option functional is activated.

See plotprocess to easily plot the coefficients and their confidence intervals.

The new algorithms have been suggested in Chernozhukov, Fernández-Val and Melly (2020a). For more detailed information about the Stata command, please refer
to Chernozhukov, Fernández-Val and Melly (2020b).

Options

Model

quantile(#) specifies the quantile(s) to be estimated.
This is either a number between 0 and 1 or a list of numbers between 0 and 1.
All the shorthands described in [U] 11.1.8 numlist can be used;
for instance quantile(0.1(0.05)0.9) is allowed.
The default value is 0,5 if the option functional is not activated.
When functional is activated, the default is a grid with meshsize 0.01.
The range of the grid depends on the sample size and the number of regressors. This avoids going too far into the tails when working with moderate sample sizes.

qlow(#),qhigh(#),qstep(#) can also be used to specify the quantiles to be estimated.
It is equivalent to quantile(qlow(qstep)qhigh). This alternative is offered only to circumvent the maximum size of a matrix or of a numeric list in Stata.
These limits, which may be an issue when a large number of quantile regressions are estimated, are not binding in Mata.

Algorithmic Method

method([mtype], [mopts])
specifies the algorithmic method used to compute the quantile regression
estimates as well as the parameters governing the algorithm.

mtype specifies the type of algorithmic method used to compute the point estimates. Available types are qreg,
fn, pqreg, pfn, proqreg, profn, and onestep.
The default method tries to select the fastest algorithm for the problem at hand.

method(qreg), selects the built-in Stata algorithm used by the official command qreg.
If several quantile regressions are being estimated, qrprocess uses the estimates of the previous quantile
regression as starting values for the next quantile regression. This already decreases the computing time compared
to sequential call of qreg.
It is based on the quantile regression estimator suggested by Koenker and Bassett (1978). The minimization problem solved by the estimator can be written as a convex linear program.
This kind of problem can be relatively efficiently solved with some form of the simplex algorithm.
The linear programming problem is a constrained optimization problem.
The feasible region is defined by the intersection of half-planes, and so it is a convex set.
The optimum for the problem must lie on the boundary of this feasible region.
The simplex begins with a feasible point at one of the vertices of the feasible set.
Then it "walks" along the edges of the feasible set from vertex to vertex, in such a way that the value of the objective function monotonically decreases at each step.
For this reason, it is an exterior point algorithm.
The simplex algorithm performs extremely well for problems of moderate size but requires as much as 150 times the computational time of OLS for a sample size of 50,000.

method(fn), selects a basic implementation of the Frisch-Newton interior point method suggested by Portnoy and Koenker (1997).
Unlike simplex, interior point algorithms actually begin in the interior of the feasible region of the linear programming program, and travels on a path towards the the boundary, converging at the optimum.
The inequality constraints are replaced by a barrier that penalizes points that are close to the boundary of the feasible set.
Since this barrier idea was pioneered by Ragnar Frisch and each iteration corresponds to a Newton step, Portnoy and Koenker (1997) call their application of the interior point method to quantile regression the Newton-Frisch algorithm.
This method is significantly faster for middle to large sample sizes.

method(pqreg) & method(pfn), selects respectively the built-in qreg exterior point method or the Frisch-Newton (fn) interior point method after the preprocessing step suggested by Portnoy and Koenker (1997).
Using preprocessing allows us to build a faster algorithm for big sample sizes.
Preprocessing uses preliminary knowledge of the sign of the residuals in order to speed up the estimation process.
The quantile regression fitted line interpolates at least k data points (where k is the number of covariates).
It can easily be shown that the quantile regression estimates are numerically identical if we change the values of the observations that are not interpolated as long as they remain on the same side of the regression line.
The sign of the residuals is the only thing that matters in the determination of the estimates.
This explains why the quantile regression estimator is influenced only by the local behavior of the conditional distribution of the response near the specified quantile and is robust to outliers in the response variable.
Following an idea in Portnoy and Koenker (1997), we exploit this property to design a quicker algorithm.
Before proceeding to the main estimations, we start by estimating an initial quantile regression using only a subsample.
This allows us to guess the sign of the residuals in the whole sample.
Using this information, we can accelerate the estimation process.

method(proqreg) and method(profn) estimates the discretized quantile regression process using respectively qreg or fn after the preprocessing step based on
the previously estimated quantile regression (preprocessing for the quantile regression process).
When many quantile regressions are estimated, then we can start with a conventional algorithm (specified by the suboption first) for the first quantile regression and then
progressively climb over the conditional quantile function, using recursively the previous quantile regression as a guess for the next one. This algorithm
provides numerically the same estimate as the traditional algorithms because we can check if the guesses were correct and stop only if this is the case.
This drastically speeds up the estimation process when the estimation is required for many quantiles.

method(1step), estimates the discretized quantile regression process using one-step iteration.
Note that, contrarily to the other algorithm, this method does not exactly minimizes the quantile objective function
of Koenker and Bassett (1978). Thus, the point estimates may differ when this method is used. However,
asymptotically, it converges to the same true if a growing number of quantile regressions is estimated.
Thus, this algorithm is intended to estimate a large number of quantile regressions.
The first regression (for the quantile that is the closest to the median) is estimated using a traditional algorithm (specified by the suboption first).
Then, sequentially, the other quantile regressions are obtained by iterating only one-time starting from the previous quantile regression.
This is the quickest method with large samples when many quantile regressions are estimated.

mopts consists of the following options: beta(#), small(#), max_it(#), m(#), and error_tol(#).
These options are technical in nature and should not be modified if the user does not understand the details of the algorithms.

beta(#) is a technical step length parameter used by the Frisch-Newton interior point method. The default value is 0.9995.

small(#) is tolerance parameter for convergence of the Frisch-Newton interior point method. The default value is 0.00001.

max_it(#) specifies the maximum number of iterations for the Frisch-Newton interior point method. The default value is 100.

m(#) is the constant m used to determine the subsample size for the methods involving preprocessing (pqreg,pfn,proqreg and profn).
The default value is 3.

error_tol(#) is an integer that determines the number of allowed mispredicted signs of residuals for the methods involving preprocessing (pqreg,pfn,proqreg and profn) along with the bootstrap.
The default value is 0.

first(algorithm) selects the algorithm for the first quantile regression when the method is proqreg, profn or 1step.
It matters only for the first quantile regression when the whole quantile regression process is estimated.
The following algorithms are available: qreg, fn, pqreg, and pfn.

Variance Covariance Matrix of the Estimates (VCE)

vce([vcetype], [vceopts])
selects the method to compute the VCE and the options related to the estimation of the variance.

vcetype selects the estimator of the VCE. Available types are kernel,nid, bootstrap, multiplier, iid and novar.

vce(kernel), selects a kernel estimate of the VCE as proposed by Powell (1991).
This is the default vcetype when functional tests and confidence bands are not requested.
This method combines relatively good performance with short computation time.
In our simulations, however, the bootstrap offered an even better performance.

vce(nid), selects the estimate of the VCE suggested by Hendricks and Koenker (1992), which use a local estimate of the sparsity.
This method is reliable but requires estimating two additional quantile regressions, and is therefore slower than kernel.

vce(bootstrap), computes bootstrap estimates of the VCE.
Estimating the VCE analytically requires the choice of a smoothing method and a bandwidth. To avoid this choice, it is popular to use resampling methods.
Chernozhukov, Fernández-Val and Melly (2013) proved the validity of the exchangeable bootstrap for the (entire) quantile regression process.
The exchangeable bootstrap covers the empirical, weighted, subsampling without replacement, and m out of n bootstraps (subsampling with replacement) as special cases.
All these special cases have been implemented in qrprocess. The specific resampling method can be determined with the suboption bmethod.
By default, the standard empirical bootstrap is used: n observations are sampled with replacement from the sample of size n.

vce(multiplier), estimates the variance by bootstrapping the first-order approximation of the quantile regression estimates.
The multiplier bootstrap takes less time to compute because it does not involve reestimating
the whole quantile regression process for each bootstrap draw, as it is the case when . The first-order approximation is computed once and is
mutiplied with independent bootstrap weights. By default, the wild bootstrap is used but Gaussian or re-centered exponential weights
can also be used as well as the empirical bootstrap or sub-sampling.

vce(iid), computes the VCE assuming iid standard errors.
We do not recommend using this method except if the researcher has very good reasons to believe in the location shift model.
In such a case all quantile regression lines are parallel. This means that this methodolgy is not consistent exactly when quantile regression is interesting.

vce(novar), does not compute the VCE. This option should be used if only the point estimates are of interest.

vceopts provides further options related to the computation of the variance.
It consists of the following options: bofinger, misspecification. cluster(varname), strata(varname), bmethod(bmethod), reps(#), subsize(#), noreplacement.

bofinger, selects Bofinger (1975) bandwidth. By default, Hall-Sheather (1988) bandwidth is used. This option is not relevant if bootstrap is chosen.

misspecification, selects the Angrist, Chernozhukov and Fernández-Val (2006) estimator of the variance, which is robust against misspecification of the conditional quantile function. This option is relevant for the kernel and nid estimators of the VCE.

cluster(varname), specifies that the standard errors allow for intragroup correlation, relaxing the usual requirement that the observations be independent. That is, the observations are independent across groups (clusters) but not necessarily within groups. The numeric cluster varname variable specifies to which group each observation belongs.

strata(varname), specifies that the sampling was stratified. This means that the population was divided into non-overlapping strata and pre-specified numbers of observations were sampled from each stratum. The numeric strata varname variable specifies to which stratum each observation belongs.

bmethod(bmethod) selects the bootstrap method. The following values for bmethod are available when vce(bootstrap) has been chosen: empirical (the default), weighted, and subsampling. The following values for bmethod are available when vce(multiplier) has been chosen: wild (the default), exponential, Gaussian, empirical, and subsampling.

empirical uses the standard empirical bootstrap. This is the default option .

subsampling uses subsampling with replacement, also called m out of n bootstrap.

weighted uses the weighted bootstrap.
For each bootstrap replication, instead of resampling the observations from the sample, all observations are kept and weighted. The weights are drawn independently from the standard exponential distribution.
One advantage of the weighted bootstrap is that no perfect mutlicolinearity problem can appear for any bootstrap replication.
Since all the bootstrap weights are always strictly positive, the rank of the regressors matrix will remain the same in each of the boostrap replication.
As a resullt, the weighted bootstrap is particularly appropriate when some regressors are indicators of rare characteristics.

wild uses wild bootstrap. This option can be used only when the multiplier bootstrap has been selected.

exponential uses recentered standard exponential weights. This option can be used only when the multiplier bootstrap has been selected.

Gaussian uses Gaussian weights. This option can be used only when the multiplier bootstrap has been selected.

reps(#), specifies the number of bootstrap replications to be performed. The default is 100.

noreplacement specifies that the resampling should be done without replacement. This method is not recommended.
It is crucial to specify a small subsample size (much smaller than the sample size) when this option is activated.
This option is relevant only if bmethod(subsampling) is chosen.

subsize(#) specifies the size of the subsamples to be drawn. This option is
relevant only if bmethod(subsampling) is chosen.

functional, activates the tools for functional inference (uniform confidence bands, test of functional hypotheses).
These tools are based on the whole quantile regression process.
Therefore, it is crucial that a large number of quantile regressions is estimated. This is done by default if the options quantile, or, qlow, qhigh, and qstep are not specified.
The tests use the Kolmogov-Smirnov and the Cramer-von-Mises statistics to measure deviations from the null hypothesis.
The p-values for five null hypotheses are computed: no effect, positive effect, negative effect, constant effect (or homogeneity hypothgesis, or location shift hypothesis) and location-scale shift.
Each of these hypothesis is tested for a single parameter or for all slope parameters simultaneously.
vce(bootstrap) or vce(multiplier) are the only accepted methods when functional is activated.

Reporting

level(#) specifies the confidence level, as a percentage, for confidence intervals. The
default is level(95) or as set by set level; see
[R] estimation options.

norpint prevents the display of the coefficients table.
By default, the coefficients table is displayed but it can be extremely voluminous when many quantile regressions are estimated.

Examples

Use the auto dataset

. sysuse auto

Median regression

. qrprocess price weight length foreign

Replay results

. qrprocess

Estimate .25 quantile regression and estimate the variance with Hendricks and Koenker (1992) method and the Bofinger bandwidth

. qrprocess price weight length foreign, quantile(.25) vce(nid, bofinger)

Median regression with bootstrap standard errors

. qrprocess price weight length foreign, vce(boot)

Estimate .25, .50, and .75 quantile regressions with bootstrap standard errors

. qrprocess price weight length foreign, vce(boot) quantile(0.25 0.5 0.75)

Estimate .9 quantile regression using subsampled bootstrap standard errors with 200 repetitions and a subsize of 50 observations

. qrprocess price weight length foreign, vce(boot, bm(subsampling) subsize(50) reps(200)) quantile(0.95)

Use the cps91 dataset

. use http://www.stata.com/data/jwooldridge/eacsap/cps91

Estimate the discretized quantile regression process for the 0.1,0.02,...,0.9 quantiles, do not print the results

. qrprocess lwage c.age##c.age i.black i.hispanic educ, quantile(0.1(0.01)0.9) noprint

Plot all the coefficients

. plotprocess

Estimate the same process, activate functional inference, use the multiplier bootstrap with 500 replications

. qrprocess lwage c.age##c.age i.black i.hispanic educ, quantile(0.1(0.01)0.9) functional vce(multiplier, reps(500))

Plot the coefficient for education with uniform and pointwise confidence bands

. plotprocess educ, ytitle("QR coefficent") title("Years of education")

Estimate the same process with the one-step estimator,

. qrprocess lwage c.age##c.age i.black i.hispanic educ, method(onestep) quantile(0.1(0.01)0.9)

Stored results

qrprocess stores the following in e():

Scalars

e(N) number of observations

e(df_m) model degrees of freedom

e(df_r) residual degrees of freedom

e(N_strat) number of strata

e(N_clust) number of clusters

e(reps) number of replications

e(subsize) subsample size

Macros

e(cmd) qrprocess

e(cmdline) command as typed

e(title) Quantile Regression

e(depvar) name of dependent variable

e(xvar) name of the regressor(s)

e(vce) vce type specified in vce()

e(bwmethod) bandwidth method; hsheather or bofinger

e(wtype) weight type

e(wexp) weight expression

e(properties) b V

e(predict) program used to implement predict

e(estat_cmd) program used to implement estat

e(method) algorithmic method

e(bmethod) bootstrap method V

e(replacement) subsampling method, "with replacement" or "without replacement"

e(clustvar) name of cluster variable

e(stratvar) name of strata variable

Matrices

let k be the number of regressors and nq be the number of quantile regressions

e(quantiles) estimated quantile(s); nq*1 vector

e(b) coefficient vector; (k*nq)*1 vector

e(V) variance-covariance matrix of the estimators; (k*nq)*(k*nq) matrix

e(sum_rdev) sum of absolute deviations; nq*1 vector

e(sum_mdev) sum of raw deviations; nq*1 vector

e(coefmat) matrix of coefficients; (k*nq) matrix

e(pointwise) lower and upper bounds of the pointwise confidence intervals; (k*nq)*2 matrix

e(uniform) lower and upper bounds of the functional confidence bands; (k*nq)*2 matrix

e(tests) p-values of the functional tests; (k + 1)*10 matrix

Functions

e(sample) marks estimation sample

Version requirements

This command requires Stata 9.2 or later. In addition, it requires the
package moremata (see Jann (2005)). Type ssc install moremata
to install this package.

References

Angrist, J., V. Chernozhukov, and I. Fernández-Val. 2006. Quantile regression under misspecification, with an application to the US wage structure. Econometrica 74: 539–563.

Bofinger, E. 1975. Estimation of a density function using order statistics. Australian Journal of Statistics 17(1): 1-7.

Chamberlain, G. 1994. Quantile regression, censoring, and the structure of wages. Advances in Econometrics, Sixth World Congress (Volume1); C. A. Sims (ed.), Cambridge University Press.

Chernozhukov, V., I. Fernández-Val, and B. Melly. 2013. Inference on counterfactual distributions. Econometrica 81(6): 2205-2268.

Chernozhukov, V., I. Fernández-Val, and B. Melly. 2020a. Fast algorithms for the quantile regression process. Working paper.

Chernozhukov, V., I. Fernández-Val, and B. Melly. 2020b. Quantile and distribution regression in Stata: algorithms, pointwise and functional inference. Working paper.

Hall, P., and S. J. Sheather. 1988. On the distribution of a studentized quantile. Journal of the Royal Statistical Society, Series B 50: 381-391.

Hendricks, W., and R. Koenker. 1992. Hierarchical spline models for conditional quantiles and the demand for electricity. Journal of the American Statistical Association 87(417): 58-68.

Jann, B. 2005. moremata: Stata module (Mata) to provide various
functions. Statistical Software Components S455001, Department of Economics,
Boston College. http://ideas.repec.org/c/boc/bocode/s455001.html.

Koenker, R. 2005. Quantile Regression. Cambridge University Press: New York.

Koenker, R., and G. Bassett. 1978. Regression quantiles. Econometrica 46: 33-50.

Portnoy, S., and R. Koenker. 1997. The Gaussian hare and the Laplacian tortoise: computability of squared-error versus absolute-error estimators. Statistical Science 12(4): 279-300.

Powell, J. L. 1991. Estimation of monotonic regression models under quantile restrictions. Nonparametric and semiparametric methods in Econometrics,(Cambridge University Press, New York) 357-384.

Remarks

This is a preliminary version. Please feel free to share your comments, reports of bugs and
propositions for extensions.

If you use this command in your work, please cite Chernozhukov, Fernández-Val and Melly (2020a) and/or
Chernozhukov, Fernández-Val and Melly (2020b).

Authors

Victor Chernozhukov, Iván Fernández-Val and Blaise Melly

MIT, Boston University and University of Bern

mellyblaise@gmail.com