[R-sig-Geo] empty neighbour sets and lm.morantest()

[Roger Bivand]

On Thu, 26 May 2011, Advait Godbole wrote:

> Dear List,
>
> I am performing a spatial analysis for gridded data with ~8000 cells. Each
> grid cell is 10 x 10 km and has one dependent and 5 independent variables.
> The dependent variable is log-transformed and has some NAs, as do the
> independent variables.
>
> I set up the lm() function using na.action=na.exclude (different from the
> default na.omit, since the lm documentation says that na.exclude pads the
> residuals to the correct length and fills with NA). Using a 1st order queen
> weights matrix generated keeping zero.policy=TRUE, I ran the lm.morantest()
> function, which returned this error:

Well, working without complete cases is always very unsatisfactory. 
Consider subsetting everything (data and neighbour object) to remove 
incomplete cases first - the results will be the same anyway. The error 
messages are caused by the internal subsetting of the neighbour object 
creating observations with no neighbours. You can choose to pass 
zero.policy=TRUE in all functions using the neighbour object, or setting 
set.ZeroPolicyOption() globally. The error message does tell you what the 
problem is.

>
>> allxform.lm.moran.queen <-
> lm.morantest(ols.allxform,vulcgrid.queen.nb.w,alternative = "two.sided")
> Error in nb2listw(neighbours = subnb, glist = NULL, style = style,
> zero.policy = zero.policy) :
>  Empty neighbour sets found
>
> repeating the same using an 8 nearest neighbour matrix also yielded the same
> error:
>
>> allxform.lm.moran.queen <-
> lm.morantest(ols.allxform,vulcgrid.knn8.w,alternative = "two.sided")
> Error in nb2listw(neighbours = subnb, glist = NULL, style = style,
> zero.policy = zero.policy) :
>  Empty neighbour sets found
>
> I then tried reading in a 45km distance thresholded neighbour matrix
> generated in OpenGeoDa using read.gwt2nb(), after converting to listw using
> zero.policy=TRUE and style="W". Running lm.morantest using this weights
> matrix returns the error:
>
>> allxform.lm.moran.queen <-
> lm.morantest(ols.allxform,vulcgrid.nb.dist.45km.gwt.w,alternative =
> "two.sided")
> Error in lm.morantest(ols.allxform, vulcgrid.nb.dist.45km.gwt.w, alternative
> = "two.sided") :
>  objects of different length
>
> Checking length of the objects, which I presume are "ols.allxform" (my lm
> result) and "vulcgrid.nb.dist.45km.gwt.w":
>
>> length(ols.allxform)
> [1] 13
>> length(vulcgrid.nb.dist.45km.gwt.w)
> [1] 3

Wrong check: length(residuals(ols.allxform)) is not equal to 
length(vulcgrid.nb.dist.45km.gwt.w$neighbours) - that is, the object 
created in OpenGeoDa is not for the same objects as the lm model objects.

Roger

>
>
> Weights matrix summary:
>
>> vulcgrid.nb.dist.45km.gwt.w
> Characteristics of weights list object:
> Neighbour list object:
> Number of regions: 8281
> Number of nonzero links: 533700
> Percentage nonzero weights: 0.7782724
> Average number of links: 64.44874
>
> Weights style: W
> Weights constants summary:
>     n       nn   S0       S1       S2
> W 8281 68574961 8281 261.9284 33173.84
>
> I then ran the lm.LMtests() function on the same linear model for the
> queen's weights matrix and for the distance thresholded weights matrix
> yileding the same error..
>
>> chi.lmtest.nohdr <- lm.LMtests(ols.allxform, vulcgrid.queen.nb.w,
> zero.policy=NULL, test="all", spChk=NULL, naSubset=TRUE)
> Error in lm.LMtests(ols.allxform, vulcgrid.queen.nb.w, zero.policy = TRUE,
> :
>  objects of different length
>
> sessionInfo output:
>
> R version 2.11.1 (2010-05-31)
> x86_64-pc-mingw32
>
> locale:
> [1] LC_COLLATE=English_United Kingdom.1252  LC_CTYPE=English_United
> Kingdom.1252
> [3] LC_MONETARY=English_United Kingdom.1252 LC_NUMERIC=C
>
> [5] LC_TIME=English_United Kingdom.1252
>
> attached base packages:
> [1] grDevices datasets  splines   graphics  stats     tcltk     utils
> methods   base
>
> other attached packages:
> [1] pgirmess_1.4.9      RandomFields_1.3.41 plyr_1.2.1
> robust_0.3-11
> [5] rrcov_1.1-00        pcaPP_1.8-3         mvtnorm_0.9-92
> robustbase_0.5-0-1
> [9] sandwich_2.2-6      zoo_1.6-4           rgdal_0.6-28
> raster_1.8-9
> [13] RNetCDF_1.2-1.1     clim.pact_2.2-41    akima_0.5-4         ncdf_1.6.3
>
> [17] car_2.0-2           nnet_7.3-1          tkrplot_0.0-19
> TeachingDemos_2.7
> [21] sphet_1.1-0         fields_6.3          spam_0.23-0         geoR_1.6-32
>
> [25] RColorBrewer_1.0-2  spdep_0.5-33        coda_0.13-5
> deldir_0.0-12
> [29] maptools_0.7-38     foreign_0.8-40      nlme_3.1-96         MASS_7.3-7
>
> [33] Matrix_0.999375-39  lattice_0.18-8      boot_1.2-42         sp_0.9-71
>
> [37] survival_2.35-8
>
> loaded via a namespace (and not attached):
> [1] grid_2.11.1   stats4_2.11.1 tools_2.11.1
>
>
> Suggestions are greatly appreciated!
>
> Regards,
>
>

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no