Jenks natural breaks classification

gistfile1.py

# code from http://danieljlewis.org/files/2010/06/Jenks.pdf
# described at http://danieljlewis.org/2010/06/07/jenks-natural-breaks-algorithm-in-python/

def getJenksBreaks( dataList, numClass ):
  dataList.sort()
  mat1 = []
  for i in range(0,len(dataList)+1):
    temp = []
    for j in range(0,numClass+1):
      temp.append(0)
    mat1.append(temp)
  mat2 = []
  for i in range(0,len(dataList)+1):
    temp = []
    for j in range(0,numClass+1):
      temp.append(0)
    mat2.append(temp)
  for i in range(1,numClass+1):
    mat1[1][i] = 1
    mat2[1][i] = 0
    for j in range(2,len(dataList)+1):
      mat2[j][i] = float('inf')
  v = 0.0
  for l in range(2,len(dataList)+1):
    s1 = 0.0
    s2 = 0.0
    w = 0.0
    for m in range(1,l+1):
      i3 = l - m + 1
      val = float(dataList[i3-1])
      s2 += val * val
      s1 += val
      w += 1
      v = s2 - (s1 * s1) / w
      i4 = i3 - 1
      if i4 != 0:
        for j in range(2,numClass+1):
          if mat2[l][j] >= (v + mat2[i4][j - 1]):
            mat1[l][j] = i3
            mat2[l][j] = v + mat2[i4][j - 1]
    mat1[l][1] = 1
    mat2[l][1] = v
  k = len(dataList)
  kclass = []
  for i in range(0,numClass+1):
    kclass.append(0)
  kclass[numClass] = float(dataList[len(dataList) - 1])
  countNum = numClass
  while countNum >= 2:#print "rank = " + str(mat1[k][countNum])
    id = int((mat1[k][countNum]) - 2)
    #print "val = " + str(dataList[id])
    kclass[countNum - 1] = dataList[id]
    k = int((mat1[k][countNum] - 1))
    countNum -= 1
  return kclass
  
def getGVF( dataList, numClass ):
  """
  The Goodness of Variance Fit (GVF) is found by taking the 
  difference between the squared deviations
  from the array mean (SDAM) and the squared deviations from the 
  class means (SDCM), and dividing by the SDAM
  """
  breaks = getJenksBreaks(dataList, numClass)
  dataList.sort()
  listMean = sum(dataList)/len(dataList)
  print listMean
  SDAM = 0.0
  for i in range(0,len(dataList)):
    sqDev = (dataList[i] - listMean)**2
    SDAM += sqDev
  SDCM = 0.0
  for i in range(0,numClass):
    if breaks[i] == 0:
      classStart = 0
    else:
      classStart = dataList.index(breaks[i])
      classStart += 1
    classEnd = dataList.index(breaks[i+1])
    classList = dataList[classStart:classEnd+1]
    classMean = sum(classList)/len(classList)
    print classMean
    preSDCM = 0.0
    for j in range(0,len(classList)):
      sqDev2 = (classList[j] - classMean)**2
      preSDCM += sqDev2
    SDCM += preSDCM
  return (SDAM - SDCM)/SDAM
  
# written by Drew
# used after running getJenksBreaks()
def classify(value, breaks):
  for i in range(1, len(breaks)):
    if value < breaks[i]:
      return i
  return len(breaks) - 1

The Jenks technique is discussed here. The blog post points to this repo for a JavaSript implementation. That repo removed their implementation of Jenks in favor of a Ckmeans.

There is an R implementation of Ckmeans in the package Ckmeans.1d.dp. See this SO question for a discussion of porting that code to Python.

Does anyone know the arguments and return values for these functions? I have the code running, but I'm not sure what the output means. For example if I run
x = [1, 1, 3, 4, 5, 5, 8, 9, 7]
getJenksBreaks(x, 3)
The output value is [0, 1, 5, 9.0]
I'm not sure how that tells me what the breaks are. The danieljlewis site forwards to http://blogs.splintdev.geog.ucl.ac.uk/, which seems to be a blank page.
Thanks!

VeloSteve, this is great that you have the code running. I still have to learn Python to do these calculations, however, I did the example you cited in RealStatistics in Excel.

The values 0, 1 and 5 that you got seems to be the lower value of the classes and 9 is the last upper value. Check the image below. Note that RealStatistics express the classes with the values existent in the table and that's why the values are not exactly the same as yours.

Best regards!

Hello everyone. I tested the code with some of my own cases and found some differences to the jenkspy library (used it previously but was looking for Jenks Breaks on pure python code recently).

Shouldn't the last while loop of the "getJenksBreaks" function include a line to make the first element of kclass equal to the first element of dataLis?

Something like:

while countNum >= 2:#print "rank = " + str(mat1[k][countNum])
    id = int((mat1[k][countNum]) - 2)
    #print "val = " + str(dataList[id])
    kclass[countNum - 1] = dataList[id]
    k = int((mat1[k][countNum] - 1))
    countNum -= 1
kclass[0] = dataList[0] #ADDED THIS LINE
return kclass

I could be wrong but this corrected any divergences between the jenkspy library and this code.

How do I use this on Google cloud

How do I use this on Google cloud