Conformal Quantile Regression

 I've written yet another vinaigrette, this time about a simple R implementation of the conformal quantile regression method introduced by Romano, Patterson and Candes  the pdf version of the vinaigrette is available here.

There, there

The LRB is distributing a daily serving of "greatest hits" from their back issues, articles that  made a big impression on readers and/or the editors.  Yesterday's article was part 1 of Derek Parfit's piece
"on the universe" which struck me as steaming pile of metaphysical nonsense.  It prompted the
following questions, with apologies to Gertrude Stein:

Why is there no there there?
Why is there there there?

Reason vs Rationality

The Gelman Blog mentioned a paper by Lorraine Daston called The Moral Economy

of Science, [Osiris, Vol. 10, Constructing Knowledge in the History of Science

(1995), pp. 2-24] which attracted my attention because long ago I read E.P.

Thompson's The Moral Economy of the English Crowd [Past & Present

No. 50 (Feb., 1971), pp. 76-136] that made an extremely persuasive case for

delicate balance between the legitimacy of the rural market for agricultural

goods and riots when this legitimacy falls apart.  Thompson prefaces his paper

with the passage from Proverbs:  "He that withholdeth Corn, the People shall

curse him: but Blessing shall be upon the Head of him that selleth it."

Daston is at pains to distinguish her thesis from Thompson's -- she is not

about the legitimacy of science in the broader culture/environment, but the

argument, which is more about the moral economy within the scientific

community, is nonetheless quite appealing.  As such departures from one's

usual reading sometimes do, this led me done a rabbit hole of related

reading concluding with a video lecture by Daston from 2010 called "Rule

of Rules:  Or how Reason became Rationality", which constitutes a cogent

historical survey of the rise of algorithms as a substitute for judgement

and reason.  Since post WWII economics has bought into the algorithmic,

mechanistic decision making zeitgeist with a vengeance, Daston's lecture

is a refreshing antidote to the familiar brainwashing of the discipline.

As I joked when I forwarded the link of the lecture to an old friend, "I'm

not sure whether the lecture is better than "Marriage Italian Style", the

last movie I've watched, but I'm working on an algorithm to decide this."

Daston has a book with Peter Galison (he of How Experiments End) on

"Objectivity", which I've not yet read, but surely must contain the following

quotation that also appears in Daston's Moral Economy paper:

    Goethe gave voice to the worries that impel mechanical objectivity when he

    preached caution in interpreting experimental results: "For here at this pass,

    this transition from empirical evidence to judgment, cognition to application,

    all the inner enemies of man lie in wait: imagination, which sweeps him away  

    on its wings before he knows his feet have left the ground; impatience; haste;

    self-satisfaction; rigidity; formalistic thought; prejudice; ease; frivolity;

    fickleness -- this whole throng and its retinue. Here they lie in ambush and

    surprise not only the active observer but also the contemplative one who

    appears safe from all passion." Johann Wolfgang Goethe, "The Experiment as

    Mediator between Object and Subject" (1792, publ. 1823)

Coincidentally, the latest Economist letters column contains related passage

from the "melancholy finale of Giuseppe Tomasi di Lampedusa’s novel, The Leopard

    Nowhere has truth such a short life as in Sicily; a fact has scarcely happened five

    minutes before its genuine kernel has vanished, been camouflaged, embellished,

    disfigured, squashed, annihilated by imagination and self interest; shame, fear,

    generosity, malice, opportunism, charity, all the passions, good as well as evil,

    fling themselves onto the fact and tear it to pieces; very soon it has vanished

    altogether

Sicily is now faced with some stiff competition in this domain.

Computational Breakthroughs for Quantile Regression

Xuming He sent me a new paper with Xiaoou Pan, Kean Ming Tan, and Wen-Xin Zhou on computation of quantile regression.  They suggest a new way to smooth the objective
function that preserves convexity and use first order coordinate descent.  At first reading I thought that their comparison might be focused exclusively on the simplex method in my quantreg package which is based on the 1974 Barrodale and Roberts algorithm, so I decided to do some comparisons myself  with the interior point method that I call the Frisch-Newton algorithm.  My initial comparisons confirm their conclusions that the new approach has massive advantages in terms of computational speed and a slight advantage in terms of accuracy.  Based on past experience I'm always skeptical about accuracy comparisons since they have often been due to some sort of shrinkage effect that while perhaps unintended had the effect of negating the affine equivariance of the original definition.  This is not the case for Xuming's new proposal.  My interpretation of the efficiency gain is that smoothing the objective function means that the new problem is in effect producing an estimator that is like local averaging of several adjacent quantiles.  One might think of it as local Huberization and it is known that provided that the bandwidth isn't too large this can yield efficiency gains.

My experiments are based on very limited experience with dense n by p designs with both X and y Gaussian.  The He et al implementation is available  in R from https://cran.r-project.org/web/packages/conquer/index.html  Sample sizes range from 1000 to 100,000, and p ranges from 10 to 200.Timings (in seconds) for the mean of 10 replications look like this:

             Frisch-Newton                   He et al

        n = 1K n = 5K n = 10K n = 100K  n = 1K n = 5K n = 10K n = 100K

p = 10   0.006  0.032   0.070    0.531  0.005  0.008   0.017    0.063

p = 50   0.051  0.292   0.499    2.056  0.010  0.016   0.030    0.228

p = 100  0.126  0.591   0.843    6.130  0.011  0.034   0.061    0.485

p = 200  0.459  1.292   2.225   21.391  0.047  0.078   0.124    0.967

Whereas mean absolute error of the coefficient estimates averaged over 10 replications is:

             Frisch-Newton                   He et al

        n = 1K n = 5K n = 10K n = 100K  n = 1K n = 5K n = 10K n = 100K

p = 10   0.438  0.149   0.103    0.034  0.414  0.146   0.096    0.033

p = 50   1.627  0.730   0.524    0.157  1.447  0.682   0.498    0.152

p = 100  3.244  1.447   1.027    0.321  2.849  1.346   0.950    0.313

p = 200  7.141  2.856   2.019    0.616  6.093  2.554   1.823    0.594

Independently, Victor Chernozhukov, Ivan Fernandez-Val and Blaise Melly have proposed new

computational methods for estimating the QR process exploiting the "globbing" ideas that Steve Portnoy and I used in our Hare and Tortoise paper some time ago.  This also yields dramatic improvements in computational speed.  For further details, see https://arxiv.org/abs/1909.05782

Taken together and combined with various new ideas for bootstrapping these new developments substantially improve the scalability of QR methods.  Congratulations, and many thanks to all

contributors to this research.

Dorfman's Group Testing

The estimable Xi'an's OG https://xianblog.wordpress.com/2020/04/09/dorfmans-group-testing/has recently pointed out that the group testing idea that has been recently put forward for Covid-19 appears as Problem IX.9.26 in Feller attributed to the late Robert Dorfman.  This problem is notable for being one of the easier problems in Feller, suitable for a serious undergraduate course and particularly apt in the present circumstances.  Dorfman's application was testing GI's for syphilis before sending them abroad into WWII.  Since I have already besmirched the reputation of Dorfman in a previous post https://davoidofmeaning.blogspot.com/search?q=dorfman it seems only fair that I recant a bit and acknowledge that this group testing idea was really quite clever and practically quite important.  I won't spoil the problem by revealing the answer, but I'll provide a snapshot of it here for those unlucky enough not to have a copy of Feller handy.  (Mine has been inherited by a terrific ex grad student who was kind enough to send me a scan.)

Empirical Bayes Confidence Intervals

The second installment of the Chamberlain (Virtual) Econometrics Seminar

featured Tim Armstrong presenting new work on Empirical Bayes Confidence

Intervals coauthored with  Michal Kolesár and Mikkel Plagborg-Møller,

and available from http://arxiv.org/abs/2004.03448.  Their approach seemed to

be motivated by recent work of Chetty using James-Stein shrinkage.  The

objective was to produce "robust" confidence intervals for such linear

shrinkage estimators, and the technique employed was a clever application

of the method of moment spaces.  This seemed to perform well in the

classical Gaussian sequence model when the mixing distribution was itself

Gaussian.  But as shown by Bruce Hansen in his discussion of the paper

following Tim's presentation, the procedure exhibits severe under

coverage in a typical Gaussian outlier model.  Thus, in the 1953 terminology

of Box it apparently possesses neither "robustness of validity", i.e. coverage,

nor "robustness of efficiency", i.e. length.

This led me to wonder how nonparametric empirical Bayes methods would

perform under such circumstances.  Bruce was kind enough to provide

code for his simulations and I spent the Easter weekend exploring this

question.  The result is yet another R Vinaigrette.  I compare two

flavors of nonparametric empirical Bayes Confidence intervals, both

based on what Efron calls G-modeling. A nonparametric estimate of the

mixing distribution called G serves as a "plug-in prior"  and intervals

are then based on posterior quantiles.  We consider two estimators for G:

the classical NPMLE of Kiefer and Wolfowitz, and the log-spline estimator

of Efron.  Spoiler alert:  in keeping with the samokritika (self-criticism)

spirit of the Vinaigrette -- Efron wins.  Links to the text and the code

for the simulations reported therein can be found at:

http://www.econ.uiuc.edu/~roger/research/ebayes/ebayes.html

~                                                                         

Ford vs. Ferrari

I recently was reminded to look again at Ignacio Sarmiento Barbieri's blog post
https://ignaciomsarmiento.github.io/2017/01/20/Quantile-and-OLS-Regression-with-Rmosek.html
that compares the speed of several implementations of quantile regression in Mosek with
that of the simplex implementation in my quantreg R package.  He finds that for problems
of modest size (n = 1000) the quantreg simplex implementation which is based on the 1972
algorithm of Barrodale and Roberts was substantial faster than any of the Mosek implementations.

However, in a more realistic econometric problem based on data from Josh Angrist (n = 97367, p = 5)
the simplex algorithm was quite a bit slower than the best of the Mosek implementations.  I wondered
how the interior point algorithms implemented in fortran in quantreg would perform.

Now, Mosek is clearly the Ferrari in this contest.  I've used Mosek in many different convex
optimization problems and always found it to be remarkably quick, responsive to special road
conditions and (especially) elegant in design.  So it requires a special sort of hubris to think
that the so-called Frisch-Newton implementations that were figuratively "made in my garage"
could compete, but lo and behold, here are the results:

Unit: milliseconds

                               expr       min        lq      mean    median

 qr.mosek.dual.rep(X, y, tau = tau)  346.9841  357.2982  365.9998  364.2814

        rq.fit.fnb(X, y, tau = tau)  193.7651  194.5455  202.1676  195.2470

       rq.fit.sfn(Xs, y, tau = tau)  330.3882  341.1280  365.0841  349.7058

         rq.fit.br(X, y, tau = tau) 7266.7237 7270.0682 7289.8437 7274.8538

        uq       max neval

  370.6428  462.1361   100

  206.7692  288.0764   100

  369.8432  451.3124   100

 7285.0870 7522.6136   100

The results for the basic (dense) Frisch-Newton appear as rq.fit.fnb while the sparse version of this
algorithm appear as rq.fit.sfn and the simplex version is rq.fit.br.