Merge branch 'develop' into update-segmentation

This commit is contained in:
Andy Eschbacher
2018-03-14 12:28:12 -04:00
189 changed files with 49373 additions and 54 deletions
+1 -1
View File
@@ -1,5 +1,5 @@
comment = 'CartoDB Spatial Analysis extension'
default_version = '0.6.1'
default_version = '0.8.1'
requires = 'plpythonu, postgis'
superuser = true
schema = cdb_crankshaft
+151 -9
View File
@@ -17,7 +17,7 @@ AS $$
num_ngbrs, permutations, geom_col, id_col)
$$ LANGUAGE plpythonu VOLATILE PARALLEL UNSAFE;
-- Moran's I Local (internal function)
-- Moran's I Local (internal function) - DEPRECATED
CREATE OR REPLACE FUNCTION
_CDB_AreasOfInterestLocal(
subquery TEXT,
@@ -27,16 +27,82 @@ CREATE OR REPLACE FUNCTION
permutations INT,
geom_col TEXT,
id_col TEXT)
RETURNS TABLE (moran NUMERIC, quads TEXT, significance NUMERIC, rowid INT, vals NUMERIC)
RETURNS TABLE (
moran NUMERIC,
quads TEXT,
significance NUMERIC,
rowid INT,
vals NUMERIC)
AS $$
from crankshaft.clustering import Moran
moran = Moran()
# TODO: use named parameters or a dictionary
return moran.local_stat(subquery, column_name, w_type,
num_ngbrs, permutations, geom_col, id_col)
result = moran.local_stat(subquery, column_name, w_type,
num_ngbrs, permutations, geom_col, id_col)
# remove spatial lag
return [(r[6], r[0], r[1], r[7], r[5]) for r in result]
$$ LANGUAGE plpythonu VOLATILE PARALLEL UNSAFE;
-- Moran's I Local (internal function)
CREATE OR REPLACE FUNCTION
_CDB_MoransILocal(
subquery TEXT,
column_name TEXT,
w_type TEXT,
num_ngbrs INT,
permutations INT,
geom_col TEXT,
id_col TEXT)
RETURNS TABLE (
quads TEXT,
significance NUMERIC,
spatial_lag NUMERIC,
spatial_lag_std NUMERIC,
orig_val NUMERIC,
orig_val_std NUMERIC,
moran_stat NUMERIC,
rowid INT)
AS $$
from crankshaft.clustering import Moran
moran = Moran()
return moran.local_stat(subquery, column_name, w_type,
num_ngbrs, permutations, geom_col, id_col)
$$ LANGUAGE plpythonu VOLATILE PARALLEL UNSAFE;
-- Moran's I Local (public-facing function)
-- Replaces CDB_AreasOfInterestLocal
CREATE OR REPLACE FUNCTION
CDB_MoransILocal(
subquery TEXT,
column_name TEXT,
w_type TEXT DEFAULT 'knn',
num_ngbrs INT DEFAULT 5,
permutations INT DEFAULT 99,
geom_col TEXT DEFAULT 'the_geom',
id_col TEXT DEFAULT 'cartodb_id')
RETURNS TABLE (
quads TEXT,
significance NUMERIC,
spatial_lag NUMERIC,
spatial_lag_std NUMERIC,
orig_val NUMERIC,
orig_val_std NUMERIC,
moran_stat NUMERIC,
rowid INT)
AS $$
SELECT
quads, significance, spatial_lag, spatial_lag_std,
orig_val, orig_val_std, moran_stat, rowid
FROM cdb_crankshaft._CDB_MoransILocal(
subquery, column_name, w_type,
num_ngbrs, permutations, geom_col, id_col);
$$ LANGUAGE SQL VOLATILE PARALLEL UNSAFE;
-- Moran's I Local (public-facing function) - DEPRECATED
CREATE OR REPLACE FUNCTION
CDB_AreasOfInterestLocal(
subquery TEXT,
@@ -132,7 +198,7 @@ AS $$
$$ LANGUAGE plpythonu VOLATILE PARALLEL UNSAFE;
-- Moran's I Local Rate (internal function)
-- Moran's I Local Rate (internal function) - DEPRECATED
CREATE OR REPLACE FUNCTION
_CDB_AreasOfInterestLocalRate(
subquery TEXT,
@@ -144,15 +210,22 @@ CREATE OR REPLACE FUNCTION
geom_col TEXT,
id_col TEXT)
RETURNS
TABLE(moran NUMERIC, quads TEXT, significance NUMERIC, rowid INT, vals NUMERIC)
TABLE(
moran NUMERIC,
quads TEXT,
significance NUMERIC,
rowid INT,
vals NUMERIC)
AS $$
from crankshaft.clustering import Moran
moran = Moran()
# TODO: use named parameters or a dictionary
return moran.local_rate_stat(subquery, numerator, denominator, w_type, num_ngbrs, permutations, geom_col, id_col)
result = moran.local_rate_stat(subquery, numerator, denominator, w_type, num_ngbrs, permutations, geom_col, id_col)
# remove spatial lag
return [(r[6], r[0], r[1], r[7], r[4]) for r in result]
$$ LANGUAGE plpythonu VOLATILE PARALLEL UNSAFE;
-- Moran's I Local Rate (public-facing function)
-- Moran's I Local Rate (public-facing function) - DEPRECATED
CREATE OR REPLACE FUNCTION
CDB_AreasOfInterestLocalRate(
subquery TEXT,
@@ -172,6 +245,75 @@ AS $$
$$ LANGUAGE SQL VOLATILE PARALLEL UNSAFE;
-- Internal function
CREATE OR REPLACE FUNCTION
_CDB_MoransILocalRate(
subquery TEXT,
numerator TEXT,
denominator TEXT,
w_type TEXT,
num_ngbrs INT,
permutations INT,
geom_col TEXT,
id_col TEXT)
RETURNS
TABLE(
quads TEXT,
significance NUMERIC,
spatial_lag NUMERIC,
spatial_lag_std NUMERIC,
orig_val NUMERIC,
orig_val_std NUMERIC,
moran_stat NUMERIC,
rowid INT)
AS $$
from crankshaft.clustering import Moran
moran = Moran()
return moran.local_rate_stat(
subquery,
numerator,
denominator,
w_type,
num_ngbrs,
permutations,
geom_col,
id_col
)
$$ LANGUAGE plpythonu VOLATILE PARALLEL UNSAFE;
-- Moran's I Rate
-- Replaces CDB_AreasOfInterestLocalRate
CREATE OR REPLACE FUNCTION
CDB_MoransILocalRate(
subquery TEXT,
numerator TEXT,
denominator TEXT,
w_type TEXT DEFAULT 'knn',
num_ngbrs INT DEFAULT 5,
permutations INT DEFAULT 99,
geom_col TEXT DEFAULT 'the_geom',
id_col TEXT DEFAULT 'cartodb_id')
RETURNS
TABLE(
quads TEXT,
significance NUMERIC,
spatial_lag NUMERIC,
spatial_lag_std NUMERIC,
orig_val NUMERIC,
orig_val_std NUMERIC,
moran_stat NUMERIC,
rowid INT)
AS $$
SELECT
quads, significance, spatial_lag, spatial_lag_std,
orig_val, orig_val_std, moran_stat, rowid
FROM cdb_crankshaft._CDB_MoransILocalRate(
subquery, numerator, denominator, w_type,
num_ngbrs, permutations, geom_col, id_col);
$$ LANGUAGE SQL VOLATILE PARALLEL UNSAFE;
-- Moran's I Local Rate only for HH and HL (public-facing function)
CREATE OR REPLACE FUNCTION
CDB_GetSpatialHotspotsRate(
+114
View File
@@ -68,6 +68,63 @@ code|quads
(52 rows)
_cdb_random_seeds
(1 row)
code|quads|diff_orig|expected|moran_stat_not_null|significance_not_null|value_comparison
01|HH|t|t|t|t|t
02|HL|t|t|t|t|t
03|LL|t|t|t|t|t
04|LL|t|t|t|t|t
05|LH|t|t|t|t|t
06|LL|t|t|t|t|t
07|HH|t|t|t|t|t
08|HH|t|t|t|t|t
09|HH|t|t|t|t|t
10|LL|t|t|t|t|t
11|LL|t|t|t|t|t
12|LL|t|t|t|t|t
13|HL|t|t|t|t|t
14|LL|t|t|t|t|t
15|LL|t|t|t|t|t
16|HH|t|t|t|t|t
17|HH|t|t|t|t|t
18|LL|t|t|t|t|t
19|HH|t|t|t|t|t
20|HH|t|t|t|t|t
21|LL|t|t|t|t|t
22|HH|t|t|t|t|t
23|LL|t|t|t|t|t
24|LL|t|t|t|t|t
25|HH|t|t|t|t|t
26|HH|t|t|t|t|t
27|LL|t|t|t|t|t
28|HH|t|t|t|t|t
29|LL|t|t|t|t|t
30|LL|t|t|t|t|t
31|HH|t|t|t|t|t
32|LL|t|t|t|t|t
33|HL|t|t|t|t|t
34|LH|t|t|t|t|t
35|LL|t|t|t|t|t
36|LL|t|t|t|t|t
37|HL|t|t|t|t|t
38|HL|t|t|t|t|t
39|HH|t|t|t|t|t
40|HH|t|t|t|t|t
41|HL|t|t|t|t|t
42|LH|t|t|t|t|t
43|LH|t|t|t|t|t
44|LL|t|t|t|t|t
45|LH|t|t|t|t|t
46|LL|t|t|t|t|t
47|LL|t|t|t|t|t
48|HH|t|t|t|t|t
49|LH|t|t|t|t|t
50|HH|t|t|t|t|t
51|LL|t|t|t|t|t
52|LL|t|t|t|t|t
(52 rows)
_cdb_random_seeds
(1 row)
code|quads
01|HH
@@ -204,6 +261,63 @@ code|quads
(52 rows)
_cdb_random_seeds
(1 row)
code|quads|diff_orig|expected|moran_stat_not_null|significance_not_null
01|HH|t|t|t|t
02|HL|t|t|t|t
03|LL|t|t|t|t
04|LL|t|t|t|t
05|LH|t|t|t|t
06|LL|t|t|t|t
07|HH|t|t|t|t
08|HH|t|t|t|t
09|HH|t|t|t|t
10|LL|t|t|t|t
11|LL|t|t|t|t
12|LL|t|t|t|t
13|HL|t|t|t|t
14|LL|t|t|t|t
15|LL|t|t|t|t
16|HH|t|t|t|t
17|HH|t|t|t|t
18|LL|t|t|t|t
19|HH|t|t|t|t
20|HH|t|t|t|t
21|LL|t|t|t|t
22|HH|t|t|t|t
23|LL|t|t|t|t
24|LL|t|t|t|t
25|HH|t|t|t|t
26|HH|t|t|t|t
27|LL|t|t|t|t
28|HH|t|t|t|t
29|LL|t|t|t|t
30|LL|t|t|t|t
31|HH|t|t|t|t
32|LL|t|t|t|t
33|HL|t|t|t|t
34|LH|t|t|t|t
35|LL|t|t|t|t
36|LL|t|t|t|t
37|HL|t|t|t|t
38|HL|t|t|t|t
39|HH|t|t|t|t
40|HH|t|t|t|t
41|HL|t|t|t|t
42|LH|t|t|t|t
43|LH|t|t|t|t
44|LL|t|t|t|t
45|LH|t|t|t|t
46|LL|t|t|t|t
47|LL|t|t|t|t
48|HH|t|t|t|t
49|LH|t|t|t|t
50|HH|t|t|t|t
51|LL|t|t|t|t
52|LL|t|t|t|t
(52 rows)
_cdb_random_seeds
(1 row)
code|quads
01|HH
+37
View File
@@ -24,6 +24,25 @@ SELECT ppoints.code, m.quads
SELECT cdb_crankshaft._cdb_random_seeds(1234);
-- Moran's I local
SELECT
ppoints.code, m.quads,
abs(avg(m.orig_val_std) OVER ()) < 1e-6 as diff_orig,
CASE WHEN m.quads = 'HL' THEN m.orig_val_std > m.spatial_lag_std
WHEN m.quads = 'HH' THEN m.orig_val_std >= 0 and m.spatial_lag_std >= 0
WHEN m.quads = 'LH' THEN m.orig_val_std < m.spatial_lag_std
WHEN m.quads = 'LL' THEN m.orig_val_std <= 0 and m.spatial_lag_std <= 0
ELSE null END as expected,
moran_stat is not null moran_stat_not_null,
significance >= 0.001 significance_not_null, -- greater than 1/1000 (default)
abs(m.orig_val - ppoints.value) <= 1e-6 as value_comparison
FROM ppoints
JOIN cdb_crankshaft.CDB_MoransILocal('SELECT * FROM ppoints', 'value') m
ON ppoints.cartodb_id = m.rowid
ORDER BY ppoints.code;
SELECT cdb_crankshaft._cdb_random_seeds(1234);
-- Spatial Hotspots
SELECT ppoints.code, m.quads
FROM ppoints
@@ -61,6 +80,24 @@ SELECT ppoints2.code, m.quads
SELECT cdb_crankshaft._cdb_random_seeds(1234);
-- Moran's I local rate
SELECT
ppoints2.code, m.quads,
abs(avg(m.orig_val_std) OVER ()) < 1e-6 as diff_orig,
CASE WHEN m.quads = 'HL' THEN m.orig_val_std > m.spatial_lag_std
WHEN m.quads = 'HH' THEN m.orig_val_std >= 0 and m.spatial_lag_std >= 0
WHEN m.quads = 'LH' THEN m.orig_val_std < m.spatial_lag_std
WHEN m.quads = 'LL' THEN m.orig_val_std <= 0 and m.spatial_lag_std <= 0
ELSE null END as expected,
moran_stat is not null moran_stat_not_null,
significance >= 0.001 significance_not_null -- greater than 1/1000 (default)
FROM ppoints2
JOIN cdb_crankshaft.CDB_MoransILocalRate('SELECT * FROM ppoints2', 'numerator', 'denominator') m
ON ppoints2.cartodb_id = m.rowid
ORDER BY ppoints2.code;
SELECT cdb_crankshaft._cdb_random_seeds(1234);
-- Spatial Hotspots (rate)
SELECT ppoints2.code, m.quads
FROM ppoints2
+160 -27
View File
@@ -1,21 +1,29 @@
"""
Moran's I geostatistics (global clustering & outliers presence)
Functionality relies on a combination of `PySAL
<http://pysal.readthedocs.io/en/latest/>`__ and the data providered provided in
the class instantiation (which defaults to PostgreSQL's plpy module's `database
access functions <https://www.postgresql.org/docs/10/static/plpython.html>`__).
"""
# TODO: Fill in local neighbors which have null/NoneType values with the
# average of the their neighborhood
import pysal as ps
from collections import OrderedDict
from crankshaft.analysis_data_provider import AnalysisDataProvider
import pysal as ps
# crankshaft module
import crankshaft.pysal_utils as pu
from crankshaft.analysis_data_provider import AnalysisDataProvider
# High level interface ---------------------------------------
class Moran(object):
"""Class for calculation of Moran's I statistics (global, local, and local
rate)
Parameters:
data_provider (:obj:`AnalysisDataProvider`): Class for fetching data. See
the `crankshaft.analysis_data_provider` module for more information.
"""
def __init__(self, data_provider=None):
if data_provider is None:
self.data_provider = AnalysisDataProvider()
@@ -28,7 +36,26 @@ class Moran(object):
Moran's I (global)
Implementation building neighbors with a PostGIS database and Moran's I
core clusters with PySAL.
Andy Eschbacher
Args:
subquery (str): Query to give access to the data needed. This query
must give access to ``attr_name``, ``geom_col``, and ``id_col``.
attr_name (str): Column name of data to analyze
w_type (str): Type of spatial weight. Must be one of `knn`
or `queen`. See `PySAL documentation
<http://pysal.readthedocs.io/en/latest/users/tutorials/weights.html>`__
for more information.
num_ngbrs (int): If using `knn` for ``w_type``, this
specifies the number of neighbors to be used to define the spatial
neighborhoods.
permutations (int): Number of permutations for performing
conditional randomization to find the p-value. Higher numbers
takes a longer time for getting results.
geom_col (str): Name of the geometry column in the dataset for
finding the spatial neighborhoods.
id_col (str): Row index for each value. Usually the database index.
"""
params = OrderedDict([("id_col", id_col),
("attr1", attr_name),
@@ -53,8 +80,38 @@ class Moran(object):
def local_stat(self, subquery, attr,
w_type, num_ngbrs, permutations, geom_col, id_col):
"""
Moran's I implementation for PL/Python
Andy Eschbacher
Moran's I (local)
Args:
subquery (str): Query to give access to the data needed. This query
must give access to ``attr_name``, ``geom_col``, and ``id_col``.
attr (str): Column name of data to analyze
w_type (str): Type of spatial weight. Must be one of `knn`
or `queen`. See `PySAL documentation
<http://pysal.readthedocs.io/en/latest/users/tutorials/weights.html>`__
for more information.
num_ngbrs (int): If using `knn` for ``w_type``, this
specifies the number of neighbors to be used to define the spatial
neighborhoods.
permutations (int): Number of permutations for performing
conditional randomization to find the p-value. Higher numbers
takes a longer time for getting results.
geom_col (str): Name of the geometry column in the dataset for
finding the spatial neighborhoods.
id_col (str): Row index for each value. Usually the database index.
Returns:
list of tuples: Where each tuple consists of the following values:
- quadrants classification (one of `HH`, `HL`, `LL`, or `LH`)
- p-value
- spatial lag
- standardized spatial lag (centered on the mean, normalized by the
standard deviation)
- original value
- standardized value
- Moran's I statistic
- original row index
"""
# geometries with attributes that are null are ignored
@@ -78,13 +135,45 @@ class Moran(object):
# find quadrants for each geometry
quads = quad_position(lisa.q)
return zip(lisa.Is, quads, lisa.p_sim, weight.id_order, lisa.y)
# calculate spatial lag
lag = ps.weights.spatial_lag.lag_spatial(weight, lisa.y)
lag_std = ps.weights.spatial_lag.lag_spatial(weight, lisa.z)
return zip(
quads,
lisa.p_sim,
lag,
lag_std,
lisa.y,
lisa.z,
lisa.Is,
weight.id_order
)
def global_rate_stat(self, subquery, numerator, denominator,
w_type, num_ngbrs, permutations, geom_col, id_col):
"""
Moran's I Rate (global)
Andy Eschbacher
Args:
subquery (str): Query to give access to the data needed. This query
must give access to ``attr_name``, ``geom_col``, and ``id_col``.
numerator (str): Column name of numerator to analyze
denominator (str): Column name of the denominator
w_type (str): Type of spatial weight. Must be one of `knn`
or `queen`. See `PySAL documentation
<http://pysal.readthedocs.io/en/latest/users/tutorials/weights.html>`__
for more information.
num_ngbrs (int): If using `knn` for ``w_type``, this
specifies the number of neighbors to be used to define the spatial
neighborhoods.
permutations (int): Number of permutations for performing
conditional randomization to find the p-value. Higher numbers
takes a longer time for getting results.
geom_col (str): Name of the geometry column in the dataset for
finding the spatial neighborhoods.
id_col (str): Row index for each value. Usually the database index.
"""
params = OrderedDict([("id_col", id_col),
("attr1", numerator),
@@ -110,8 +199,39 @@ class Moran(object):
def local_rate_stat(self, subquery, numerator, denominator,
w_type, num_ngbrs, permutations, geom_col, id_col):
"""
Moran's I Local Rate
Andy Eschbacher
Moran's I Local Rate
Args:
subquery (str): Query to give access to the data needed. This query
must give access to ``attr_name``, ``geom_col``, and ``id_col``.
numerator (str): Column name of numerator to analyze
denominator (str): Column name of the denominator
w_type (str): Type of spatial weight. Must be one of `knn`
or `queen`. See `PySAL documentation
<http://pysal.readthedocs.io/en/latest/users/tutorials/weights.html>`__
for more information.
num_ngbrs (int): If using `knn` for ``w_type``, this
specifies the number of neighbors to be used to define the spatial
neighborhoods.
permutations (int): Number of permutations for performing
conditional randomization to find the p-value. Higher numbers
takes a longer time for getting results.
geom_col (str): Name of the geometry column in the dataset for
finding the spatial neighborhoods.
id_col (str): Row index for each value. Usually the database index.
Returns:
list of tuples: Where each tuple consists of the following values:
- quadrants classification (one of `HH`, `HL`, `LL`, or `LH`)
- p-value
- spatial lag
- standardized spatial lag (centered on the mean, normalized by the
standard deviation)
- original value (roughly numerator divided by denominator)
- standardized value
- Moran's I statistic
- original row index
"""
# geometries with values that are null are ignored
# resulting in a collection of not as near neighbors
@@ -138,7 +258,20 @@ class Moran(object):
# find quadrants for each geometry
quads = quad_position(lisa.q)
return zip(lisa.Is, quads, lisa.p_sim, weight.id_order, lisa.y)
# spatial lag
lag = ps.weights.spatial_lag.lag_spatial(weight, lisa.y)
lag_std = ps.weights.spatial_lag.lag_spatial(weight, lisa.z)
return zip(
quads,
lisa.p_sim,
lag,
lag_std,
lisa.y,
lisa.z,
lisa.Is,
weight.id_order
)
def local_bivariate_stat(self, subquery, attr1, attr2,
permutations, geom_col, id_col,
@@ -177,12 +310,12 @@ class Moran(object):
def map_quads(coord):
"""
Map a quadrant number to Moran's I designation
HH=1, LH=2, LL=3, HL=4
Input:
@param coord (int): quadrant of a specific measurement
Output:
classification (one of 'HH', 'LH', 'LL', or 'HL')
Map a quadrant number to Moran's I designation
HH=1, LH=2, LL=3, HL=4
Args:
coord (int): quadrant of a specific measurement
Returns:
classification (one of 'HH', 'LH', 'LL', or 'HL')
"""
if coord == 1:
return 'HH'
@@ -192,17 +325,17 @@ def map_quads(coord):
return 'LL'
elif coord == 4:
return 'HL'
else:
return None
return None
def quad_position(quads):
"""
Produce Moran's I classification based of n
Input:
@param quads ndarray: an array of quads classified by
1-4 (PySAL default)
Output:
@param list: an array of quads classied by 'HH', 'LL', etc.
Map all quads
Args:
quads (:obj:`numpy.ndarray`): an array of quads classified by
1-4 (PySAL default)
Returns:
list: an array of quads classied by 'HH', 'LL', etc.
"""
return [map_quads(q) for q in quads]
+1 -1
View File
@@ -1,5 +1,5 @@
joblib==0.8.3
numpy==1.6.1
scipy==0.14.0
pysal==1.11.2
pysal==1.14.3
scikit-learn==0.14.1
+1 -1
View File
@@ -41,7 +41,7 @@ setup(
# The choice of component versions is dictated by what's
# provisioned in the production servers.
# IMPORTANT NOTE: please don't change this line. Instead issue a ticket to systems for evaluation.
install_requires=['joblib==0.8.3', 'numpy==1.6.1', 'scipy==0.14.0', 'pysal==1.11.2', 'scikit-learn==0.14.1'],
install_requires=['joblib==0.8.3', 'numpy==1.6.1', 'scipy==0.14.0', 'pysal==1.14.3', 'scikit-learn==0.14.1'],
requires=['pysal', 'numpy', 'sklearn'],
@@ -71,10 +71,10 @@ class MoranTest(unittest.TestCase):
random_seeds.set_random_seeds(1234)
result = moran.local_stat('subquery', 'value',
'knn', 5, 99, 'the_geom', 'cartodb_id')
result = [(row[0], row[1]) for row in result]
result = [(row[0], row[6]) for row in result]
zipped_values = zip(result, self.moran_data)
for ([res_val, res_quad], [exp_val, exp_quad]) in zipped_values:
for ([res_quad, res_val], [exp_val, exp_quad]) in zipped_values:
self.assertAlmostEqual(res_val, exp_val)
self.assertEqual(res_quad, exp_quad)
@@ -89,11 +89,11 @@ class MoranTest(unittest.TestCase):
moran = Moran(FakeDataProvider(data))
result = moran.local_rate_stat('subquery', 'numerator', 'denominator',
'knn', 5, 99, 'the_geom', 'cartodb_id')
result = [(row[0], row[1]) for row in result]
result = [(row[0], row[6]) for row in result]
zipped_values = zip(result, self.moran_data)
for ([res_val, res_quad], [exp_val, exp_quad]) in zipped_values:
for ([res_quad, res_val], [exp_val, exp_quad]) in zipped_values:
self.assertAlmostEqual(res_val, exp_val)
def test_moran(self):