# Why are the predictions of a quantile regression model changed by an increasing transformation of the DV?

Take a dataset and suppose we fit two quantile regression models to it, one with the raw dependent variable (DV) and one with the logged DV. Then look at each model's predictions for the training data, reversing the transformation in the case of the logged model. Here's an example in R:

library(quantreg)

set.seed(1)

x1 = rnorm(100)
y = exp(x1 + 3*rnorm(100))

m = rq(y ~ x1, tau = .5)
p1 = predict(m)

m = rq(log(y) ~ x1, tau = .5)
p2 = exp(predict(m))



The two models give different predictions:

        1         2         3         4         5         6
0.4631787 0.8619277 0.3602179 1.5567692 0.9337257 0.3676802
1         2         3         4         5         6
0.2776789 0.7859069 0.2122640 4.8162168 0.9478212 0.2164372


But how can this be? The models are fit in terms of the quantiles of the DV, and the natural logarithm is an increasing function, so it preserves quantile ranks.

The equivariance to monotone transformations property that implies $Q_q(y \vert x)=\exp \{Q_q(\ln y|x )\}$ is exact only if the conditional quantile function is correctly specified. This is unlikely to be the case in practice, and is not the case in your simulation, since $\exp \{x+\varepsilon\} \ne x + \varepsilon$. The only case where the linear model will be exact is when all regressors are discrete and we specify a fully saturated model with dummy variables as regressors that exhaust all the possible interactions.

This will give you much better results, since I am adding a constant to avoid undefined logs of zeros and negatives, instead of exponentiating:

library(quantreg)
set.seed(1)
x1 = rnorm(100)
y = 10 + x1 + 3*rnorm(100)
m = rq(y ~1 + x1, tau = .5)
p1 = predict(m)
m2 = rq(log(y) ~1 + x1, tau = .5)
p2 = exp(predict(m2))

• What exactly is $Q_q$ here? – Kodiologist Aug 8 '16 at 19:22
• It's the $q$th conditional quantile of $y$ given $x$. The conditional median from your model would be $p1=Q_{0.5}(y \vert x1)$. It is analogous to conditional mean that you get from regression, $E[y \vert x]$. – Dimitriy V. Masterov Aug 8 '16 at 19:28

This isn't a complete answer, but it may help.

The covariate in your model may be a red herring, the issue is reproducible with an intercept only model

set.seed(1)

y = runif(100, 0, 1)  # Need to keep it positive for the log.

m = rq(y ~ 1, tau = .5)
p1 = predict(m)

m2 = rq(log(y) ~ 1, tau = .5)
p2 = exp(predict(m2))



which results in

> print(head(p1))
1         2         3         4         5         6
0.4820801 0.4820801 0.4820801 0.4820801 0.4820801 0.4820801
1         2         3         4         5         6
0.4935413 0.4935413 0.4935413 0.4935413 0.4935413 0.4935413


But if we redo this whole thing with an odd number of data points

set.seed(1)

y = runif(101, 0, 1)  # <- Now it's odd!

m = rq(y ~ 1, tau = .5)
p1 = predict(m)

m2 = rq(log(y) ~ 1, tau = .5)
p2 = exp(predict(m2))



we get

    1         2         3         4         5         6
0.4935413 0.4935413 0.4935413 0.4935413 0.4935413 0.4935413
1         2         3         4         5         6
0.4935413 0.4935413 0.4935413 0.4935413 0.4935413 0.4935413


Which is on the nose.

It's easy to reproduce this without the quantile regression call

> x <- c(1, 2, 3, 4, 5, 6)
> median(x)
 3.5
> exp(median(log(x)))
 3.464102
>
> x <- c(1, 2, 3, 4, 5, 6, 7)
> median(x)
 4
> exp(median(log(x)))
 4


I suspect this has something to do with the behaviour you're seeing, but it's not clear how to incorporate the covariate into the argument.

• Seems like a complete answer to me: because n is even and the median is computed as the average of two transformed values. – xan Aug 8 '16 at 19:57
• That is, (f(a)+f(b))/2 != f((a+b)/2) – xan Aug 8 '16 at 20:25

The fitted values are calculated by substituting in the vector values of the independent variable to the coefficient estimates (in this case of the median).

$y = mx +c$.

There is a direct linear relationship between the original (m) and log (m2) predictions (ie. same x with different slope and intercept).

Here's a code snippet to get to the nuts and bolts of it:

set.seed(1)
x1 = rnorm(100)
y = exp(x1 + 3*rnorm(100))
X <- as.matrix(cbind("(Intercept)"=rep(1,length(x1)),x1))
rq1 <- rq.fit.fnb(X, y, tau = 0.5)
resid <- (y - X %*% rq1$coefficients) fit <- X %*% rq1$coefficients

rq2 <- rq.fit.fnb(X, log(y), tau = 0.5)
residlog <- (y - X %*% rq2$coefficients) fitlog <- X %*% rq2$coefficients

model <- coef(lm(fit~fitlog))
plot(fitlog ~ fit)
[,1]
[1,] 0.4769196
[2,] 0.8292802
[3,] 0.3859367
[4,] 1.4432873
[5,] 0.8927255
[6,] 0.3925308

[,1]
[1,] 0.4769196
[2,] 0.8292802
[3,] 0.3859367
[4,] 1.4432873
[5,] 0.8927255
[6,] 0.3925308


Refs for the mathematics of quantile regression: