Friday, September 16, 2016

Dawn of the δ-method





Several years ago my colleague Steve Portnoy wrote a letter to the editor of the American Statistician
in response to an article that they had published called Who invented the δ-method?  The article
claimed priority for Robert Dorfman on the basis of an article appearing in 1938 called "A Note
on the δ-method for Finding Variance Formulae" published in the Biometric Bulletin.  Portnoy
pointed out that Joe Doob had written about the δ-method in a 1935 Annals paper titled,  "On the
limiting distribution of certain statistics" referring to it as the "well-known δ-method" and citing
prior work by T.L. Kelley and Sewell Wright, and noting rather modestly that his Theorem 1
"shows an interpretation which can be given to the results obtained by this method."  It seems
plausible that Doob's is the first formal justification for the method, and it is puzzling to put it
euphemistically  that Dorfman made no mention of Doob's article.  Perhaps this oversight can be
forgiven as a juvenile mistake since the Dorfman paper was written shortly after he finished his
undergraduate studies at Columbia, while working at the Worcester State Hospital, pictured above.
This august institution was reputed to be the first asylum for the insane in New England, and also happened to be the publisher of the Biometric Bulletin.  Dorfman later went on to earn a Phd at Berkeley, and taught at Harvard where hecoauthored an influential book about linear programming with Paul Samuelson and Robert Solow.

Thursday, September 15, 2016

Bag of Little Bootstrap for QR

In the never ending quest to speed up inference for large quantile regression problems,
I have started to look into the Kleiner, et al Bag of Little Bootstraps proposal.  After a
serious confusion on my part was corrected with the help of Xiaofeng Shao,  I've come
to the following code fragment added to summary.rq in my quantreg package:


 else if (se == "BLB"){ # Bag of Little Bootstraps
        n <- length(y)
        b <- ceiling(n^gamma)
        S <- n %/% b
        U <- matrix(sample(1:n, b * S), b, S)
        Z <- matrix(0, NCOL(x), S)
        for(i in 1:S){
            u <- U[,i]
           B <- matrix(0, NCOL(x), R)
           for(j in 1:R){
w <- c(rmultinom(1, n, rep(1/b, b)))
B[,j] <- rq.wfit(x[u,], y[u], tau, weights = w, method = "fnb")$coef
           }
            Z[,i] <- sqrt(diag(cov(t(B))))
        }
        serr <- apply(Z, 1, mean)
    }

In the eventual implementation I managed to embed the inner loop into fortran
which helps to speed things up a bit, although of course it would be eventually
helpful to allow this to be distributed across cluster nodes.
I should also mention that the implicit assumption here that BLB works for
moments, appears to be (presently) beyond the scope of the current theory,