Roger,

Out of curiosity, did you ever reach a final word on this?

Adam

On Fri, Feb 24, 2012 at 7:18 AM, Roger Bivand <Roger.Bivand@nhh.no> wrote:

> On Fri, 24 Feb 2012, Adam Slez wrote:
>
> > Am I correct in thinking that the degrees of freedom associated with the
> > bptest.sarlm, much like the bptest code from which it came, should be
> > driven by the number of parameters in the auxiliary regression NOT the
> > model used to generate the residuals?  The way the code for bptest.sarlm
> is
> > written, the degrees of freedom are calculated from the original
> > regression, as opposed to the auxiliary model.  The easiest way to see
> this
> > is to do something like this:
> >
> > library(spdep)
> > data(columbus)
> > wgts<-nb2listw(col.gal.nb)
> > model1<-lagsarlm(CRIME~INC, data = columbus, listw = wgts)
> > model2<-lagsarlm(CRIME~INC+HOVAL, data = columbus, listw = wgts)
> > bp1<-bptest.sarlm(model1, ~EW, data = columbus, studentize = F)
> > bp2<-bptest.sarlm(model2, ~EW, data = columbus, studentize = F)
> > bp3<-bptest.sarlm(model2, ~EW+CP+NSA, data = columbus, studentize = F)
> >
> > All I'm doing here is running spatially adjusted BP tests on two
> different
> > spatial lag models.  In the first two cases, I'm just testing for
> groupwise
> > heteroskedasticity with regimes defined in terms of the EW variable.  If
> I
> > understand the logic of the BP test, the degrees of freedom associated
> with
> > bp1 and bp2 should be the same, while the number of degrees of freedom
> > associated with bp2 and bp3 should differ.  This doesn't seem to be the
> > case.
> >
> > If we repeat this exercise for non-spatial models and tests, changes in
> the
> > degrees of freedom work as expected (recognizing of course that the
> actual
> > degrees of freedom for a non-spatial BP tests differs from the degrees of
> > freedom for a spatially-adjusted BP test):
> >
> > library(lmtest)
> > model3<-lm(CRIME~INC, data = columbus)
> > model4<-lm(CRIME~INC+HOVAL, data = columbus)
> > bp4<-bptest(model3, ~EW, data = columbus, studentize = F)
> > bp5<-bptest(model4, ~EW, data = columbus, studentize = F)
> > bp6<-bptest(model4, ~EW+CP+NSA, data = columbus, studentize = F)
> >
> > The bptest.sarlm routine is cribbed directly from the bptest code.  If we
> > take a look a the first part of the bptest.sarlm code, we can see the
> > problem (I'm using ... to denote skipped code):
> >
> > function (object, varformula = NULL, studentize = TRUE, data = list()) {
> >    if (!inherits(object, "sarlm"))
> >        stop("not sarlm object")
> >    Z <- object$tarX
> >    k <- ncol(Z)
> >    n <- nrow(Z)
> >    if (!is.null(varformula))
> >        Z <- model.matrix(varformula, data = data)
>
> Thanks, k and n should be set after the if(), not before, I think. I'll
> check and release if this seems correct.
>
> Roger
>
> >   ...
> >    df <- k - 1
> >   ...
> > }
> >
> > It is pretty clear that df is driven by k which is unaffected by Z <-
> > model.matrix(varformula, data = data), even when varformula is not null.
>  I
> > can't see any reason why the degrees of freedom for a spatially-adjusted
> BP
> > test would follow from the regression producing the residuals as opposed
> to
> > the auxiliary regression.  Any help would be much appreciated.
> >
> > Adam
> >
> >       [[alternative HTML version deleted]]
> >
> > _______________________________________________
> > R-sig-Geo mailing list
> > R-sig-Geo@r-project.org
> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> >
>
> --
> Roger Bivand
> Department of Economics, NHH Norwegian School of Economics,
> Helleveien 30, N-5045 Bergen, Norway.
> voice: +47 55 95 93 55; fax +47 55 95 95 43
> e-mail: Roger.Bivand@nhh.no
>
>

	[[alternative HTML version deleted]]

