Optim() by default is using Nelder-Mead which is an extremely poor way to
do linear programming, despite the fact that ?optim says that: "It will work
reasonably well for
non-differentiable functions." I didn't check your coding of the objective
function fully, but at the
very least you should explicitly pass the arguments y, x, and quant. and you
need to replace
what you call sample.cond.quantile by 0 in the definition of w.less.
url: www.econ.uiuc.edu/~roger Roger Koenker
email rkoen...@uiuc.edu Department of Economics
vox: 217-333-4558 University of Illinois
fax: 217-244-6678 Urbana, IL 61801
On Jun 7, 2012, at 1:49 PM, Kevin Chang wrote:
> Hello Everyone,
>
>
>
> I'm currently learning about quantile regressions. I've been using an
> optimizer to compare with the rq() command for quantile regression.
>
> When I run the code, the results show that my coefficients are consistent
> with rq(), but the intercept term can vary by a lot.
>
> I don't think my optimizer code is wrong and suspects it has something to do
> with the starting values.
>
> The results seems very sensitive to different starting values and I don't
> know how to make sense of it.
>
>
>
> Advice from the community would be greatly appreciated.
>
>
>
> Sincerely,
>
>
> Kevin Chang
>
>
>
> ###################### CODE Below ###########################
>
>
>
> library(quantreg)
>
> data(engel)
>
> y<-cbind(engel[,2])
>
> x<-cbind(rep(1,length(engel[,1])),engel[,1])
>
> x1<-cbind(engel[,1])
>
> nn<-nrow(engel)
>
> nn
>
>
>
> bhat.ls<-solve(t(x)%*%x)%*%t(x)%*%y
>
> #bhat.ls
>
>
>
> # QUANTILES
>
> quant=.25
>
>
>
> fr.1=function(bhat.opt)
>
> {
>
> uu=y-x%*%bhat.opt
>
> sample.cond.quantile=quantile(uu,quant)
>
> w.less=rep(0,nn)
>
> for(ii in 1:nn){if(uu[ii]<sample.cond.quantile) w.less[ii]=1}
>
>
>
> sum((quant-1)*sum((y-x%*%bhat.opt)*w.less) #negative residuals
>
> +quant*sum((y-x%*%bhat.opt)*(1-w.less))) #positive residuals
>
> }
>
> start<-c(0,0)
>
> result=optim(start,fr.1)
>
> bhat.cond=result$par
>
>
>
> #Quantile Command Results
>
> fit.temp=rq(y~x1,tau=quant)
>
> fit.temp
>
>
>
> #OPTIMIZER Results
>
> bhat.cond
>
>
>
> #OLS Command Results
>
> mean=lm(y~x1)
>
> mean
>
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help@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@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.
