Spatial Autocorrelation with rasters

[HolmSeifert]

Hallo,

Can anyone tell me a procedure to get a spatial autocorrelation based on gridcells.
Here is the problem.
I want to use the morans autocorrelation for a raster dataset.
The formula or equation is clear for me.
But the core problem is now to calculate the cross product for a defined neighborhood.
The goal is to calculate local morans autocorrelation values for every grid cell. I guess, this is not easy as it sounds.
My consideration was to calculate all the values from the focal statistics. The variance and summed weighted values, which are needed in the morans calculation are no problems.
But the cross product in the numerator make some difficulties.
I found no solution to get the cross product for a 3x3 neighborhood or higher.

Can anyone help me.

Holm

[JeffreyEvans]

The attached zip file contains an aml that will return the Moran's-I coefficient using either a Rook or Queen's case. It runs in Workstation ArcInfo Grid. I cannot perform this analysis in ArcMap because there is nothing equivalent to a DOCELL in Spatial Analyst. I also have code that will calculate Geary's-C.

[PhilMorefield]

It runs in Workstation ArcInfo Grid. I cannot perform this analysis in ArcMap because there is nothing equivalent to a DOCELL in Spatial Analyst.

For the most part, the equivalent to DOCELL is Map Algebra/Raster Calculator. It's the ability to use neighborhood notation that is missing.

With ArcGIS 10, you could convert your raster to a Numpy array and do the calculations that way.

Also, you could just convert your raster to points or polygons, couldn't you?

[HolmSeifert]

Hallo,

First thanks to all the people who wrote a reply.

1) Jeffrey Evans,
thanks for your aml file. Unfortunately i do no have a arcinfo workstation.
So i have to perform / create this analysis in arcgis itself.

2) Phil Morefield,
Sure this would be a solution to run a morans analysis in python.
Converting to points is no problems. But my task is to solve such a operation with grid cells.
The reasons are due to many thematic maps existing in raster format.

3) William Huber
Thanks for your very detailled solution.
I have some question.
But before i start i will provide a simple example.

Example.
Let´s imagine a simple 3x3 grid with a indicator transformation to 1 or 0.
I had done this before with my rasters. This looks like

1 1 0
0 1 0
0 0 1

For a better understanding i will give names to every cell.

a b c
d e f
g h i

My weight w(ij) looks like this:
0.5  0.5  0.5
0.5  1.0  0.5
0.5  0.5  0.5

Ok no i want to calculate the morans

So every gridcell should have its own morans.
"I" stands for morans between -1 and +1
The result should be like this

Ia Ib Ic
Id Ie If
Ig Ih Ii

I have studied your equation.

I have attached an example.
It is an gdb geodatabase and a weighted kernel.

Big letters explains the source or intermediate results and small letters are the operators.
TestGrid is the base.



Ok lets beginn. Here are the grids.


I Denumerator
"Denumerator1_Variance" = is the focal (3x3) variance of the values
"Denumerator2_Variance_times_DoubleSumWeights" = is the product of "Denumerator1_Variance" and the Doublesum of the 3x3 weight. The doublesum of the weight is 10 (weight matrix above)

So the denumerator is solved.

But now the numerator.


II Numerator.
"Numerator1_focalXj" = is the Sum over all neighboring gridcells from exactly one center cell with the weights / weightmatrix.


"Numerator2_FocalXj_minus_X = is the subtraction of the value normal x-value from the "Numerator1_FocalXj" grid.


"Numerator3_Numerator2_product_TestGrid = is the product of the "Numerator2_FocalXj_minus_X" grid and the normal "TestGrid"


"Numerator4_focal_Numerator3" = is the FocalSum of the "Numerator3_product_Numerator2" grid.
!!!I used here a focal sum Willima instead of a zonal sum, because i want the morans for every grid cell in a 3x3 neighborhood. i think this should be the correct procedure.


"Numerator5_Numerator4_times_8" = is the product of the "Numerator4_focal_Numerator3" grid and the number of neighbors, here 8 neighbors.

So we have the numerator.
http://forums.arcgis.com/images/icons/icon6.png



III Divide
The last calculation is: Numerator / Denumerator

"Divide" = "Numerator4_focal_Numerator3" grid divide by the "Denumerator2_Variance_times_DoubleSumWeights" grid.

this should be the endresult. For me it looks ok.

Please test and write me back.

Thanks a lot

Holm

[PhilMorefield]

As I mentioned previously ArcGIS is capable of cell-by-cell calculations, but the "neighborhood notation" available in ArcInfo no longer exists.

I haven't followed the necessary calculations in detail, but it seems like some combination of Block Statistics and Focal Statistics should also work.

Edit:

Also, in ArcGIS 10 neighborhood clases (e.g. rectange, wedge, circle) can be defined explicitly.