adds hotspot/coldspot function

This commit is contained in:
Andy Eschbacher 2016-09-09 11:11:32 -04:00
parent 2261d11de0
commit 60f52633fa
3 changed files with 68 additions and 0 deletions

16
src/pg/sql/16_getis.sql Normal file
View File

@ -0,0 +1,16 @@
-- Getis-Ord's G
-- Hotspot/Coldspot Analysis tool
CREATE OR REPLACE FUNCTION
CDB_GetisOrdsG(
subquery TEXT,
column_name TEXT,
w_type TEXT,
num_ngbrs INT,
permutations INT,
geom_col TEXT,
id_col TEXT)
RETURNS TABLE (z_val NUMERIC, p_val NUMERIC, rowid BIGINT)
AS $$
from crankshaft.clustering import getis
return getis_ord(subquery, column_name, w_type, num_ngbrs, permutations, geom_col, id_col)
$$ LANGUAGE plpythonu;

View File

@ -1,3 +1,4 @@
"""Import all functions from for clustering"""
from moran import *
from kmeans import *
from getis import *

View File

@ -0,0 +1,51 @@
"""
Moran's I geostatistics (global clustering & outliers presence)
"""
# TODO: Fill in local neighbors which have null/NoneType values with the
# average of the their neighborhood
import pysal as ps
import plpy
from collections import OrderedDict
# crankshaft module
import crankshaft.pysal_utils as pu
# High level interface ---------------------------------------
def getis_ord(subquery, attr,
w_type, num_ngbrs, permutations, geom_col, id_col):
"""
Getis-Ord's G
Implementation building neighbors with a PostGIS database and Getis-Ord's G
hotspot/coldspot analysis with PySAL.
Andy Eschbacher
"""
# geometries with attributes that are null are ignored
# resulting in a collection of not as near neighbors if kNN is chosen
qvals = OrderedDict([("id_col", id_col),
("attr1", attr),
("geom_col", geom_col),
("subquery", subquery),
("num_ngbrs", num_ngbrs)])
query = pu.construct_neighbor_query(w_type, qvals)
try:
result = plpy.execute(query)
# if there are no neighbors, exit
if len(result) == 0:
return pu.empty_zipped_array(3)
except plpy.SPIError, e:
plpy.error('Query failed: %s' % e)
attr_vals = pu.get_attributes(result)
weight = pu.get_weight(result, w_type, num_ngbrs)
# calculate LISA values
getis = ps.esda.getisord(attr_vals, weight, star=True)
return zip(getis.z_sim, getis.p_sim, weight.id_order)