R-squared in quantile regression

Question

I am using quantile regression to find predictors of 90th percentile of my data. I am doing this in R using the quantreg package. How can I determine $r^2$ for quantile regression which will indicate how much of variability is being explained by predictor variables?

What I really want to know: "Any method I can use to find how much of variability is being explained?". Significance levels by P values is available in output of command: summary(rq(formula,tau,data)). How can I get goodness of fit?

Answer 1

Koenker and Machado$^{[1]}$ describe $R^1$, a local measure of goodness of fit at the particular ($\tau$) quantile.

Let $V(\tau) = \min_{b}\sum \rho_\tau(y_i-x_i'b)$

Let $\hat{\beta}(\tau)$ and $\tilde{\beta}(\tau)$ be the coefficient estimates for the full model, and a restricted model, and let $\hat{V}$ and $\tilde{V}$ be the corresponding $V$ terms.

They define the goodness of fit criterion $R^1(\tau) = 1-\frac{\hat{V}}{\tilde{V} }$.

Koenker gives code for $V$ here,

rho <- function(u,tau=.5)u*(tau - (u < 0))
V <- sum(rho(f$resid, f$tau))

So if we compute $V$ for a model with an intercept-only ($\tilde{V}$ - or V0 in the code snippet below) and then an unrestricted model ($\hat{V}$), we can calculate an R1 <- 1-Vhat/V0 that's - at least notionally - somewhat like the usual $R^2$.

Edit: In your case, of course, the second argument, which would be put in where f$tau is in the call in the second line of code, will be whichever value of tau you used. The value in the first line merely sets the default.

'Explaining variance about the mean' is really not what you're doing with quantile regression, so you shouldn't expect to have a really equivalent measure.

I don't think the concept of $R^2$ translates well to quantile regression. You can define various more-or-less analogous quantities, as here, but no matter what you choose, you won't have most of the properties real $R^2$ has in OLS regression. You need to be clear about what properties you need and what you don't -- in some cases it may be possible to have a measure that does what you want.

--

$[1]$ Koenker, R and Machado, J (1999),
Goodness of Fit and Related Inference Processes for Quantile Regression,
Journal of the American Statistical Association, 94:448, 1296-1310

Answer 2

The pseudo-$R^2$ measure suggested by Koenker and Machado (1999) in JASA measures goodness of fit by comparing the sum of weighted deviations for the model of interest with the same sum from a model in which only the intercept appears. It is calculated as

$$R_1(\tau) = 1 - \frac{\sum_{y_i \ge \hat y_i} \tau \cdot \vert y_i-\hat y_i \vert +\sum_{y_i<\hat y_i} (1-\tau) \cdot \vert y_i-\hat y_i \vert}{\sum_{y_i \ge \bar y} \tau \cdot \vert y_i-\bar y \vert +\sum_{y_i<\bar y_i} (1-\tau) \cdot \vert y_i-\bar y \vert},$$

where $\hat y_i =\alpha_{\tau}+\beta_{\tau}x$ is the fitted $\tau$th quantile for observation $i$, and $\bar y=\beta_{\tau}$ is the fitted value from the intercept-only model.

$R_1(\tau)$ should lie in $[0,1]$, where 1 would correspond to a perfect fit since the numerator which consists of the weighted sum of deviations would be zero. It a local measure of fit for QRM since it depends on $\tau$, unlike the global $R^2$ from OLS. That is arguably the source of the warnings about using it: if you model fits in the tail, there's not guarantee that it fits well anywhere else. This approach could also be used to compare nested models.

Here's an example in R:

library(quantreg)
data(engel)

fit0 <- rq(foodexp~1,tau=0.9,data=engel)
fit1 <- rq(foodexp~income,tau=0.9,data=engel)

rho <- function(u,tau=.5)u*(tau - (u < 0))
R1 <- 1 - fit1$rho/fit0$rho

This could probably be accomplished more elegantly.