[R-sig-Geo] How to handle NA-values in raster-based Geary´s C test?

Kerstin Traut Kerstin.Traut at uni-jena.de
Thu Jul 12 19:31:04 CEST 2012 

Zitat von Roger Bivand <Roger.Bivand at nhh.no>:

> On Thu, 12 Jul 2012, Kerstin Traut wrote:
>
>> Zitat von Roger Bivand <Roger.Bivand at nhh.no>:
>>
>>> On Thu, 12 Jul 2012, Kerstin Traut wrote:
>>>
>>>> Zitat von Roger Bivand <Roger.Bivand at nhh.no>:
>>>>
>>>>> On Thu, 12 Jul 2012, Kerstin Traut wrote:
>>>>>
>>>>>> Zitat von Roger Bivand <Roger.Bivand at nhh.no>:
>>>>>>
>>>>>>> On Thu, 12 Jul 2012, Roger Bivand wrote:
>>>>>>>
>>>>>>>> On Thu, 12 Jul 2012, Kerstin Traut wrote:
>>>>>>>>
>>>>>>>>> Zitat von Roger Bivand <Roger.Bivand at nhh.no>:
>>>>>>>>>
>>>>>>>>>> On Thu, 12 Jul 2012, Kerstin Traut wrote:
>>>>>>>>>>
>>>>>>>>>>> Hi,
>>>>>>>>>>>
>>>>>>>>>>> I have a question on testing spatial autocorrelation on  
>>>>>>>>>>> raster data including NA-values. In particular, I like to  
>>>>>>>>>>> calculate Moran´s I and Geary´s C indices by using inverse  
>>>>>>>>>>> distance weighting matrices.
>>>>>>>>>>> Calculating Moran´s I with moran.test works fine, because  
>>>>>>>>>>> the function contains the option "na.action=na.pass".  
>>>>>>>>>>> Unfortunately, the function geary.test does not contain  
>>>>>>>>>>> this option. The problem is that geary.test needs an nb  
>>>>>>>>>>> argument for which I cannot perform for example  
>>>>>>>>>>> nb<-ifelse(value, nb, NA), because it would change the  
>>>>>>>>>>> nb-format to a list-format.
>>>>>>>>>>> Is there a way to ignore NA-values in geary.test? Or how  
>>>>>>>>>>> can I set NA values in a neigbors (nb) list?
>>>>>>>>>>
>>>>>>>>>> The simple answer is to subset the variable of interest and  
>>>>>>>>>> the nb object to omit NA observations (which is what  
>>>>>>>>>> moran.test does internally). Alternatively, you could  
>>>>>>>>>> contribute a patch to geary.test(), if you like.
>>>>>>>>>>
>>>>>>>>>> Hope this helps,
>>>>>>>>>>
>>>>>>>>>> Roger
>>>>>>>>>
>>>>>>>>> Thanks for your answer. Unfortunately, I think that  
>>>>>>>>> subsetting is in my case not the right solution, because I  
>>>>>>>>> need to define the amount of neighbors in cell2nb on the  
>>>>>>>>> basis of a rectangular extent (e.g. 50 x 50 neighbors).  
>>>>>>>>> Subsetting this extent by deleting NA-values would destroy  
>>>>>>>>> the nb list format required by geary.test.
>>>>>>>>> Is there any solution for subsetting an nb list without  
>>>>>>>>> destroying the nb format?
>>>>>>>>
>>>>>>>> Yes, there are subset methods for nb and listw objects. To  
>>>>>>>> ensure data integrity, use dnearneigh() instead of cell2nb(),  
>>>>>>>> and use a SpatialPointsDataFrame to hold your data. In that  
>>>>>>>> way, you make sure that the subset condition is applied to  
>>>>>>>> the same observations.
>>>>>>>
>>>>>>> The latter advice may not work if your coordinates are planar  
>>>>>>> and you really mean to set torus=TRUE in cell2nb(). You'll  
>>>>>>> have to get the data order right by deconstructing the  
>>>>>>> region.id attribute of the nb object created by cell2nb() to  
>>>>>>> ensure data integrity.
>>>>>>>
>>>>>>> Roger
>>>>>>
>>>>>> I already tested dnearneigh(), but I like to avoid this  
>>>>>> function because for large datasets cell2nb() calculates  
>>>>>> neighbors more quickly.
>>>>>> Thanks to your suggestion I was able to subset the nb list by  
>>>>>> using nb.subset(). But I could not find a suitable function for  
>>>>>> subsetting a listw. Using the subset() function results in the  
>>>>>> following error message: "Not yet able to subset general  
>>>>>> weights lists". Do you have an explaination?
>>>>>
>>>>> Yes, because the changes may not be appropriate, depending on  
>>>>> how the general weights were defined. In your case, subset the  
>>>>> nb object and the list of idw weights, and rebuild the listw  
>>>>> object with the subsetted objects.
>>>>>
>>>>> Roger
>>>>
>>>> That´s what I already tried, but the error message says  
>>>> "neighbors and glist do not conform"??
>>>
>>> Right, you'll need to rebuild from the subsetted nb object from  
>>> the nbdists() step (and subsetted coordinates), because the  
>>> subsetted idw list will have weights for excluded observations.
>>>
>>> Roger
>>
>> Ok, I hope I understood everything right. You mean, firstly, I need  
>> to subset my nb list and afterwards, I need to calculate the  
>> nb-distances using nbdists()?
>
> Yes, that will give consistent results when the listw object  
> contains general weights. The same would apply to moran.test() with  
> a listw object with general weights, I would have thought based on  
> the code of subset.listw(), so I suggest not using na.pass in the  
> moran.test() call, and doing the same reconstruction before calling  
> the tests. I'd be interested in seeing an example of data (value,  
> nb.idw, coordinates(subset)),  that run successfully through  
> moran.test() because it should also fail on general weights.
>
> Roger

Finally, after some more failed trials the geary.test() is working for  
my dataset. Thank you so much for your suggestions!
I tested moran.test() with the consistent method against the na.pass  
method. For the first method I got an I-value of 0.7952, for the  
second a Moran´s I value of 0.79477. Do you know, how to explain this  
small (but existing) difference?

Kerstin

>>
>> Kerstin
>>
>>
>>>>
>>>> Kerstin
>>>>
>>>>>>
>>>>>> Kerstin
>>>>>>
>>>>>>
>>>>>>>>
>>>>>>>> Roger
>>>>>>>>
>>>>>>>>>
>>>>>>>>> Kind regards,
>>>>>>>>> Kerstin
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>>> Thank you very much,
>>>>>>>>>>> Kerstin
>>>>>>>>>>>
>>>>>>>>>>> ---
>>>>>>>>>>>
>>>>>>>>>>> file<-raster(test.tif")
>>>>>>>>>>> mask<-extent(file, 1, 50, 1, 50)
>>>>>>>>>>> subset<-crop(file,mask,test.tif", overwrite=TRUE)
>>>>>>>>>>> value <- as.matrix(subset)
>>>>>>>>>>>
>>>>>>>>>>> # Creation of a list of integer vectors giving the region  
>>>>>>>>>>> id numbers for the neighbors within the grid extent
>>>>>>>>>>> nb <- cell2nb(50, 50, type="queen", torus=TRUE)
>>>>>>>>>>>
>>>>>>>>>>> # Calculation of the distances along the links in the  
>>>>>>>>>>> neighbous list
>>>>>>>>>>> dist <- nbdists(nb, coordinates(subset), longlat=TRUE)
>>>>>>>>>>>
>>>>>>>>>>> # Converting the distances to inverse distance values
>>>>>>>>>>> inv.dist <- lapply(dist, function (x) 1/(x*100))
>>>>>>>>>>>
>>>>>>>>>>> # Supplementing the neighbors list with spatial weights  
>>>>>>>>>>> for a coding scheme
>>>>>>>>>>> nb.idw <- nb2listw(nb, glist=inv.dist, style="W")
>>>>>>>>>>> summary(unlist(nb.idw$weights))
>>>>>>>>>>>
>>>>>>>>>>> # Moran's test (two-sided)
>>>>>>>>>>> moran.global <- moran.test(value, listw=nb.idw,  
>>>>>>>>>>> zero.policy = TRUE, na.action=na.pass,  
>>>>>>>>>>> alternative="two.sided")
>>>>>>>>>>> moran.global
>>>>>>>>>>>
>>>>>>>>>>> # Geary's test (two-sided)
>>>>>>>>>>> geary.global <- geary.test(value, listw=nb.idw,  
>>>>>>>>>>> zero.policy = TRUE , alternative="two.sided")
>>>>>>>>>>> geary.global
>>>>>>>>>>>
>>>>>>>>>>> ----------------------------------------------------------------
>>>>>>>>>>> This message was sent through https://webmail.uni-jena.de
>>>>>>>>>>>
>>>>>>>>>>> _______________________________________________
>>>>>>>>>>> R-sig-Geo mailing list
>>>>>>>>>>> R-sig-Geo at 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 at nhh.no
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> ----------------------------------------------------------------
>>>>>>>>> This message was sent through https://webmail.uni-jena.de
>>>>>>>>
>>>>>>>>
>>>>>>>
>>>>>>> -- 
>>>>>>> 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 at nhh.no
>>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> ----------------------------------------------------------------
>>>>>> This message was sent through https://webmail.uni-jena.de
>>>>>
>>>>> -- 
>>>>> 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 at nhh.no
>>>>>
>>>>
>>>>
>>>>
>>>> ----------------------------------------------------------------
>>>> This message was sent through https://webmail.uni-jena.de
>>>
>>> -- 
>>> 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 at nhh.no
>>>
>>
>>
>>
>> ----------------------------------------------------------------
>> This message was sent through https://webmail.uni-jena.de
>
> -- 
> 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 at nhh.no
>



----------------------------------------------------------------
This message was sent through https://webmail.uni-jena.de

More information about the R-sig-Geo mailing list