Quantile regression and p-values - getting more decimal places

Question

Using R, and package quantreg, I am performing quantile regression analyses to my data.

I can get access to the p-values using the se (standard error) estimator in the summary function, as below, however I only get 5 decimal places, and would like more.

model <- rq(outcome ~ predictor)
summary(model, se="ker")

Call: rq(formula = outcome ~ predictor)

tau: [1] 0.5

Coefficients:
            Value    Std. Error t value  Pr(>|t|)
(Intercept) 78.68182  2.89984   27.13312  0.00000
predictor    0.22727  0.03885    5.84943  0.00000

How might I get access to more decimal places on the p-values?

---

Update

Ok, so I can get some more decimal places by selecting the sub-object that contains the matrix of numerical results;

> summary(model, se="ker")[[3]]
                 Value Std. Error   t value     Pr(>|t|)
(Intercept) 78.6818182 3.13897835 25.066059 0.000000e+00
predictor    0.2272727 0.04105681  5.535567 4.397638e-08

However the P-value is still rounded to 0 when the value is <1e-12 (the above output is a simplified example model). I can get some more by applying the suggestion from @seancarmody ;

format(summary(model, se="ker")[[3]], digits=22)

But if P < 1e-22 it is still rounded to 0, and if "digits" is set to > 22 I get the following error;

format(summary(model, se="ker")[[3]], digits=23)

Error in prettyNum(.Internal(format(x, trim, digits, nsmall, width, 3L, : invalid 'digits' argument

Is it possible to access even more decimal places?

Answer 1

To get any farther I think you have to dig in and see how the p values are calculated. In particular, summary.rq has the following snippet:

  coef[, 4] <- if (rdf > 0) 
        2 * (1 - pt(abs(coef[, 3]), rdf))
    else NA

This is actually a fairly imprecise calculation of the p-value (which under ordinary circumstances doesn't really matter). You can probably get the maximum amount of precision by retrieving the log of the p-value [for example, you could in principle retrieve p-values less than 10^{-308}, the smallest value that R can represent as a double-precision value], e.g.

ss <- summary(model,se="ker")
log(2)+pt(abs(ss$coefficients[,"t value"]),
     lower.tail=FALSE,log.p=TRUE,df=ss$rdf)

The lower.tail=FALSE argument gives you the complement (upper-tail) value of the CDF; log.p=TRUE says you want the log value; adding the log(2) makes it two-sided.

Answer 2

Have a look at str(model). You can see there is an attribute coefficients, which will give you a numeric vector with higher precision than displayed by summary. You can also look at these with

format(model$coefficients, digits=12)

Note that this converts the coefficients to strings.

To examine the p.values more closely, try

summary(model, se="ker")$coefficients[,1]*1e27

although I very much doubt that there is any meaning to be found in the digits of a p.value less than 1e-27!