R for Spatial Statistics

 

Cluster Analysis

Just in case you wish to try a Mantel test here is an R function that calculates a permutated Mantel and associated p-value. You do not need any additional packages. The input are  two matrix objects. There is an option to plot the permutation.  Just copy/paste the function into R and call it using mantel. There is an example in the header. I imagine that it would be fairly easy to translate this code to Python so it could be called via ArcGIS but there would need to be some coercion of the data into distance matrices.  
 
##########################################################################
# PROGRAM: MantelPermutation. R
# USE: MANTEL TEST OF TWO MATRICES VIA PERMUTATION
# REQUIRES: NA
#           
# ARGUMENTS: 
#       q1      MATRIX 1
#       q2      MATRIX 2
#       nperm   NUMBER OF PERMUTATIONS
#       plot    OPTION TO PLOT PERMUTATION
#       ...     ARGUMENTS PASSED TO PLOT 
#
# VALUE:
#       A LIST OBJECT WITH MANTEL z AND p VALUE
#
# EXAMPLES: 
#      ( q1 = matrix(runif(36),nrow=6) )
#      ( q2 = matrix(runif(36),nrow=6) )
#        mantel(q1, q2, plot=TRUE, nperm=199, rpt=200) 
#
# CONTACT: 
#     Jeffrey S. Evans 
#     Senior Landscape Ecologist 
#     The Nature Conservancy - Central Science
#     Adjunct Assistant Professor
#     Zoology and Physiology
#     University of Wyoming
#     Laramie, WY
#     (970)672-6766
#     jeffrey_evans@tnc.org
##########################################################################
mantel = function(m1, m2, nperm=250, plot=TRUE, rpt=FALSE) {
    factorial = function(n)gamma(n+1)
    perm.rowscols = function(m1,n,p=NULL) {
       if (is.numeric(p))
         s = gen.perm(p,n)
       else
         s = sample(1:n)
       m1[s,s]
    }   
      mant.zstat = function(m1, m2) { sum(lower.triang(m1*m2)) }   
    gen.perm = function(k,n,vec=1:n) {
       p = as.list(vec)
       perm = numeric(n)
       for (i in (n:1)) {
          v = (k%%i)+1
          perm[i] = p[[v]]
          p[[v]] = NULL
          k = k %/% i
       }
      perm
    }  
    all.perm = function(n,vec=1:n) {
       k = matrix(NA,nrow=factorial(n),ncol=n)
       for (i in (1:factorial(n)))
           k[i,] = gen.perm(i,n,vec)
       k
    }   
    lower.triang = function(m) {
       d = dim(m)
       if (d[1] != d[2]) print("Warning: non-square matrix")
       m[col(m)factorial(n)) {
      cat("Warning: nperm>n!, enumerating permutations =",factorial(n),"\n")
      nperm = factorial(n)
      enum=T
   }
   nullstats = numeric(nperm)
   for (i in (1:nperm)) {
      if (enum)
         nullstats[i] = mant.zstat(m1,perm.rowscols(m2,n,i))
      else
         nullstats[i] = mant.zstat(m1,perm.rowscols(m2,n))
      if (is.numeric(rpt))
         if (i %% rpt == 0)
            cat(i,"...\n")
   }
   nullstats = sort(nullstats)
   pval = 1-((1:nperm)[nullstats>realz])[1]/nperm
   if (plot==TRUE) {
      plot(density(nullstats),type="l",
       main="Distribution of Mantel z-statistic",
       xlab="Z statistic",ylab="# of permutations",
       sub=paste("Actual z-stat =",round(realz,3),
         ": p<",round(pval,4),",",nperm," permutations"))
      abline(v=realz)
   }
   return( list(Z.statistic=realz,pValue=pval) )
}
######################## END FUNCTION ################################
 
Jeffrey S. Evans
Senior Landscape Ecologist  
The Nature Conservancy
Central Science/DbyD
Adjunct Assistant Professor
Zoology & Physiology 
University of Wyoming
Laramie, WY 82070 
jeffrey_evans@tnc.org
(970) 672-6766