I'm working on a data set consisting of data about geographical districts. Some of the variables are continuous for which I manage to calculate global and local autocorrelation using R package spdep. The data set also has key variables in 9-level ordinal scale. I'm a confused about how to calculate autocorrelation for such data. I have understood that at least for binary nominal data (such as presence/absence) joint count -statistic would be a fitting method, but can it be extented to ordinal or multilevel nominal data?

Any help is greatly appreciated!