roger koenker wrote:>
> John,
>
> You can make a local version of bootcov which either:
>
> deletes these arguments from the call to fitter, or
> modify the switch statement to include rq.fit,
>
> the latter would need to also modify rq() to return a fitFunction
> component, so the first option is simpler. One of these days I'll
> incorporate clustered se's into summary.rq, but meanwhile
> this seems to be a good alternative.
>
> Roger
With some luck I'll have a Design package interface to quantreg for the
case where one value of tau is requested, along with a corresponding
addition to bootcov. With some luck I'll have this by tomorrow. This
would give you nomograms, effect plots, LaTeX expression of model fits, etc.
As a side issue, as I'm just getting into quantile regression, I've done
some simulations to better understand the operating characteristics.
I've seen that the efficiency of rq estimates at fixed covariate values,
when you have a normal distribution, relative to least squares
estimates, is about the same as the efficiency of the sample median
compared to the sample mean, or a bit worse. Code is below. As is
always the case, if you are willing to assume a lot (e.g., normality
with constant variance) you can do better. On the other hand, quantile
regression is more invariant (but not perfectly so) to transformations
of Y as studied using the 2nd bit of code below.
On another issue, I had trouble installing quantreg on one Ubuntu
machine, as I've forgotten which LAPACK/ATLAS ubuntu/debian package
needs to be installed. I got an error related to the absence of
llapack. If anyone has a pointer I'd appreciate getting it.
This gives me an opportunity to thank Roger for his wonderful work in
this area, and his excellent vignette that goes with quantreg.
Frank Harrell
n <- 50
m <- 1000
a <- b <- aq <- bq <- numeric(m)
for(i in 1:m)
{
cat(i,'\r')
x1 <- rnorm(n)
y <- x1 + rnorm(n)/4
f <- lm(y ~ x1)
k <- coef(f)
a[i] <- k[1]; b[i] <- k[2]
f <- rq(y ~ x1)
k <- coef(f)
aq[i] <- k[1]; bq[i] <- k[2]
}
cat('\n')
z <- c(mean(a^2), mean(aq^2))
zq <- c(mean((b-1)^2), mean((bq-1)^2))
c(z, z[1]/z[2])
c(zq, zq[1]/zq[2])
-------------------
n <- 50
par(mfrow=c(5,4), mar=c(2,1,1,1))
for(i in 1:20)
{
x1 <- rnorm(n)
y <- exp(x1 + rnorm(n)/4)
range(y)
f <- lm(y ~ exp(x1))
x1s <- seq(-2.25, 2.25, length=100)
d <- data.frame(x1=x1s)
plot(x1s, log(predict(f,d)), type='n', ylab='',
xlim=c(-2.5,2.5), ylim=c(-4, 3.5))
lines(x1s, log(predict(f,d)), col='blue')
k <- coef(f)
a <- k[1]; b <- k[2]
c(a^2, 2*a*b, b^2)
f <- lm(log(y) ~ x1)
lines(x1s, predict(f,d), col='green')
f <- rq(y ~ exp(x1))
# plot(x1s, log(predict(f,d)), type='n', ylab='')
lines(x1s, log(predict(f,d)), col='pink')
k <- coef(f)
a <- k[1]; b <- k[2]
c(a^2, 2*a*b, b^2)
f <- rq(log(y) ~ x1)
lines(x1s, predict(f,d), col='red')
}
>
> url: www.econ.uiuc.edu/~roger Roger Koenker
> email rkoenker at uiuc.edu Department of Economics
> vox: 217-333-4558 University of Illinois
> fax: 217-244-6678 Urbana, IL 61801
>
>
>
> On Jul 23, 2009, at 8:10 PM, John Gardner wrote:
>
>> I have a quick question, and I apologize in advance if, in asking, I
>> expose my woeful ignorance of R and its packages. I am trying to use
>> the bootcov function to estimate the standard errors for some
>> regression quantiles using a cluster bootstrap. However, it seems that
>> bootcov passes arguments that rq.fit doesn't like, preventing the
>> command from executing. Here is an example:
>>
>> e<-bootcov(rq(y~x),clust,B=10,fitter=rq.fit)
>>
>> (where clust is my clustering variable) results in
>>
>> Error in rq.fit.br(x, y, tau = tau, ...) :
>> unused argument(s) (maxit = 15, penalty.matrix = NULL)
>>
>> In contrast, the lm.fit function seems to just ignore these arguments,
>> resulting in the following warning:
>>
>> 10: In fitter(X[obs, , drop = FALSE], Y[obs, , drop = FALSE], maxit =
>> maxit, :
>> extra arguments maxitpenalty.matrix are just disregarded.
>>
>> Is there a way that I can either (a) modify bootcov so that it
doesn't
>> pass these arguments or (b) modify rq so that it ignores them?
>>
>> Thanks in advance,
>> John Gardner
>>
>> ______________________________________________
>> R-help at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
--
Frank E Harrell Jr Professor and Chair School of Medicine
Department of Biostatistics Vanderbilt University