# [R] ANOVA in quantreg - faulty test for 'nesting'?

Roger Koenker rkoenker at illinois.edu
Sat Apr 21 17:00:31 CEST 2012

```Yes, the models are nested, and yes I probably should have been more clever
about parsing formulae, for these cases.  I'll have a look at the code for glm, I
presume that lm() is also capable of groking this?   Meanwhile, I would like to
point out that any of these linear restrictions can be reformulated as an exclusion
restriction, see for example the remark on pages 5-6 of

http://www.econ.uiuc.edu/~roger/research/ranks/ranktests.pdf

or the references cited there.  And using this sort of parameterization you can
use anova.rq().

Roger Koenker
rkoenker at illinois.edu

On Apr 20, 2012, at 7:59 PM, Galen Sher wrote:

> Thanks to Brian's suggestion of using the logLik() function, I've dug a
> little deeper. I definitely think f1 and f2 are nested models. For example,
> by adding x2 to fmla1, we obtain a formula that quite clearly nests fmla1
> and achieves the same log likelihood as that obtained for f2. Here is the
> extra code to show this:
>
> fmla3 = y~I(x1+x2)+x2
> f3=glm(fmla3)
> logLik(f1); logLik(f2); logLik(f3)
>
> If f2=f3, as the log likelihoods would suggest, then this gives us a
> workaround: define the intermediate formula fmla3 and the fit f3 as above,
> and then conduct the analysis of variance between models f1 and f3 instead
> of f1 and f2. This doesn't offend anova.rq() any longer:
>
> f3.qr = rq(fmla3)
> anova(f1.qr,f3.qr) #this is actually anova(f1.qr, f2.qr) which resulted in
> an error earlier
>
> -Galen
>
> On Fri, Apr 20, 2012 at 6:47 PM, Brian S Cade <cadeb at usgs.gov> wrote:
>
>> Galen:  Interesting, first time I tried it I couldn't get anova.glm to
>> compute the p-values (got a warning) but I tried it again and it worked.
>> Your larger (alternative hypothesis) model is  y = B0 + B1*x1 + B2*x2 + e
>> and your smaller (null hypotheisis) model is y = B0 + B3*(x1 + x2).  So  I
>> guess I see that what you're trying to test in this case is that B1 = B2.
>> I don't think Roger Koenker anticipated such a test with anova.rq.  Other
>> options besides using information criteria (AIC, BIC, etc) for comparing
>> nonnested models include the Vuong test.  But not sure how readily the
>> theory of Vuong's test (like a paired t-test) extends to quantile
>> regression.
>>
>> Brian
>>
>>
>> U. S. Geological Survey
>> Fort Collins Science Center
>> 2150 Centre Ave., Bldg. C
>> Fort Collins, CO  80526-8818
>>
>> tel:  970 226-9326
>>
>>
>> From: Galen Sher <galensher at gmail.com> To: Brian S Cade <cadeb at usgs.gov>
>> Date: 04/20/2012 09:59 AM Subject: Re: [R] ANOVA in quantreg - faulty
>> test for 'nesting'?
>> ------------------------------
>>
>>
>>
>> Thanks Brian. I think anova.glm() requires the user to specify the
>> appropriate distribution. In the example above, if I use either of the
>> following commands
>>
>> anova(f1,f2,test="Chisq") #or
>> anova(f1,f2,test="F")
>>
>> then anova.glm() will compute and display p-values associated with the
>> deviance statistics. The reason I thought these models are nested is
>> because the first model can be thought of as the second model estimated
>> under an additional linear equality constraint. I suppose that's less of an
>> R question and more of a regression question.
>>
>> Thanks for the logLik suggestion. In the absence of more information I'll
>> have to do that - I'm just wary of conducting the test myself!
>>
>> Regards,
>> Galen
>>
>> On Fri, Apr 20, 2012 at 4:31 PM, Brian S Cade <*cadeb at usgs.gov*<cadeb at usgs.gov>>
>> wrote:
>> Galen:  Actually don't see how the models are nested (ask yourself what
>> parameter in the model with 3 parameters is set to zero in the 2 parameter
>> model?) and indeed if I try your code anova.glm will compute the difference
>> in deviance between the two models but it does not compute a probability
>> value for that difference in deviance as that would not make sense for
>> models that aren't nested.   Koenker's implementation of anova.rq
>> immediately detects that the models aren't nested so doesn't even compute
>> the deviance difference.  You could use the logLik function on the rq
>> objects to get their log likelihoods or use AIC (BIC) to compare the
>> quantile regression models.
>>
>> Brian
>>
>>
>> U. S. Geological Survey
>> Fort Collins Science Center
>> 2150 Centre Ave., Bldg. C
>> Fort Collins, CO  80526-8818
>>
>> tel:  *970 226-9326* <970%20226-9326>
>>
>>  From: galen <*galensher at gmail.com* <galensher at gmail.com>>  To: *
>> r-help at r-project.org* <r-help at r-project.org>  Date: 04/19/2012 02:40 PM
>> Subject: [R] ANOVA in quantreg - faulty test for 'nesting'?  Sent by: *
>> r-help-bounces at r-project.org* <r-help-bounces at r-project.org>
>>
>> ------------------------------
>>
>>
>>
>>
>> I am trying to implement an ANOVA on a pair of quantile regression models
>> in
>> R. The anova.rq() function performs a basic check to see whether the models
>> are nested, but I think this check is failing in my case. I think my models
>> are nested despite the anova.rqlist() function saying otherwise. Here is an
>> example where the GLM ANOVA regards the models as nested, but the quantile
>> regression ANOVA tells me the models aren't nested:
>>
>> y = rnorm(100)
>> x1 = rnorm(100)
>> x2 = rnorm(100)
>>
>> fmla1 = y~I(x1+x2)
>> fmla2 = y~x1+x2
>>
>> f1 = glm(fmla1)
>> f2 = glm(fmla2)
>>
>> anova(f1,f2) #This works
>>
>> f1.qr = rq(fmla1)
>> f2.qr = rq(fmla2)
>>
>> anova(f1.qr,f2.qr) #Error!
>> #Error in anova.rqlist(object, ...) : Models aren't nested
>>
>> Are the models in fact not nested? If they are nested, is there an easy
>> workaround I could use? Many thanks in anticipation.
>>
>> --
>> View this message in context: *
>> http://r.789695.n4.nabble.com/ANOVA-in-quantreg-faulty-test-for-nesting-tp4571994p4571994.html
>> *<http://r.789695.n4.nabble.com/ANOVA-in-quantreg-faulty-test-for-nesting-tp4571994p4571994.html>
>> Sent from the R help mailing list archive at Nabble.com.
>>
>> ______________________________________________*
>> **R-help at r-project.org* <R-help at r-project.org> mailing list*
>> **https://stat.ethz.ch/mailman/listinfo/r-help*<https://stat.ethz.ch/mailman/listinfo/r-help>
>> http://www.R-project.org/posting-guide.html*<http://www.r-project.org/posting-guide.html>
>> and provide commented, minimal, self-contained, reproducible code.
>>
>>
>>
>>
>>
>
> 	[[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help