cartodb4/crankshaft
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Merge branch 'develop' into spatial_lag

Browse Source
This commit is contained in:
Andy Eschbacher
2018-03-02 09:00:40 -05:00
parent 09255b586b b547d35882
commit 39bb6884a3
73 changed files with 16298 additions and 55 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
src/pg/crankshaft.control
View File

@@ -1,5 +1,5 @@
comment = 'CartoDB Spatial Analysis extension'
default_version = '0.6.1'
default_version = '0.7.0'
requires = 'plpythonu, postgis'
superuser = true
schema = cdb_crankshaft

59
src/pg/sql/11_kmeans.sql
View File

@@ -1,18 +1,58 @@
-- Spatial k-means clustering
CREATE OR REPLACE FUNCTION CDB_KMeans(query text, no_clusters integer, no_init integer default 20)
RETURNS table (cartodb_id integer, cluster_no integer) as $$
CREATE OR REPLACE FUNCTION CDB_KMeans(
query TEXT,
no_clusters INTEGER,
no_init INTEGER DEFAULT 20
)
RETURNS TABLE(
cartodb_id INTEGER,
cluster_no INTEGER
) AS $$
from crankshaft.clustering import Kmeans
kmeans = Kmeans()
return kmeans.spatial(query, no_clusters, no_init)
from crankshaft.clustering import Kmeans
kmeans = Kmeans()
return kmeans.spatial(query, no_clusters, no_init)
$$ LANGUAGE plpythonu VOLATILE PARALLEL UNSAFE;
-- Non-spatial k-means clustering
-- query: sql query to retrieve all the needed data
-- colnames: text array of column names for doing the clustering analysis
-- no_clusters: number of requested clusters
-- standardize: whether to scale variables to a mean of zero and a standard
-- deviation of 1
-- id_colname: name of the id column
CREATE OR REPLACE FUNCTION CDB_WeightedMeanS(state Numeric[],the_geom GEOMETRY(Point, 4326), weight NUMERIC)
RETURNS Numeric[] AS
$$
CREATE OR REPLACE FUNCTION CDB_KMeansNonspatial(
query TEXT,
colnames TEXT[],
no_clusters INTEGER,
standardize BOOLEAN DEFAULT true,
id_col TEXT DEFAULT 'cartodb_id'
)
RETURNS TABLE(
cluster_label text,
cluster_center json,
silhouettes numeric,
inertia numeric,
rowid bigint
) AS $$
from crankshaft.clustering import Kmeans
kmeans = Kmeans()
return kmeans.nonspatial(query, colnames, no_clusters,
standardize=standardize,
id_col=id_col)
$$ LANGUAGE plpythonu VOLATILE PARALLEL UNSAFE;
CREATE OR REPLACE FUNCTION CDB_WeightedMeanS(
state NUMERIC[],
the_geom GEOMETRY(Point, 4326),
weight NUMERIC
)
RETURNS Numeric[] AS $$
DECLARE
newX NUMERIC;
newY NUMERIC;
@@ -32,7 +72,8 @@ BEGIN
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION CDB_WeightedMeanF(state Numeric[])
CREATE OR REPLACE FUNCTION CDB_WeightedMeanF(state NUMERIC[])
RETURNS GEOMETRY AS
$$
BEGIN

44
src/pg/sql/13_PIA.sql
View File

@@ -31,7 +31,7 @@ DECLARE
sqr numeric;
p geometry;
BEGIN
sqr := |/2;
sqr := 0.5*(|/2.0);
polygon := ST_Transform(polygon, 3857);
-- grid #0 cell size
@@ -46,6 +46,7 @@ BEGIN
SELECT array_agg(c) INTO cells FROM c1;
-- 1st guess: centroid
best_c := polygon;
best_d := cdb_crankshaft._Signed_Dist(polygon, ST_Centroid(Polygon));
-- looping the loop
@@ -56,6 +57,7 @@ BEGIN
EXIT WHEN i > n;
cell := cells[i];
i := i+1;
-- cell side size, it's square
@@ -63,13 +65,14 @@ BEGIN
-- check distance
test_d := cdb_crankshaft._Signed_Dist(polygon, ST_Centroid(cell));
IF test_d > best_d THEN
best_d := test_d;
best_c := cells[i];
best_c := cell;
END IF;
-- longest distance within the cell
test_mx := test_d + (test_h/2 * sqr);
test_mx := test_d + (test_h * sqr);
-- if the cell has no chance to contains the desired point, continue
CONTINUE WHEN test_mx - best_d <= tolerance;
@@ -94,29 +97,46 @@ END;
$$ language plpgsql IMMUTABLE PARALLEL SAFE;
-- signed distance point to polygon with holes
-- negative is the point is out the polygon
-- rev 1. adding MULTIPOLYGON and GEOMETRYCOLLECTION support by @abelvm
CREATE OR REPLACE FUNCTION _Signed_Dist(
IN polygon geometry,
IN point geometry
)
RETURNS numeric AS $$
DECLARE
pols geometry[];
pol geometry;
i integer;
j integer;
within integer;
w integer;
holes integer;
dist numeric;
d numeric;
BEGIN
dist := 1e999;
SELECT LEAST(dist, ST_distance(point, ST_ExteriorRing(polygon))::numeric) INTO dist;
SELECT CASE WHEN ST_Within(point,polygon) THEN 1 ELSE -1 END INTO within;
SELECT ST_NumInteriorRings(polygon) INTO holes;
IF holes > 0 THEN
FOR i IN 1..holes
LOOP
SELECT LEAST(dist, ST_distance(point, ST_InteriorRingN(polygon, i))::numeric) INTO dist;
END LOOP;
END IF;
WITH collection as (SELECT (ST_dump(polygon)).geom as geom) SELECT array_agg(geom) into pols FROM collection;
FOR j in 1..array_length(pols, 1)
LOOP
pol := pols[j];
d := dist;
SELECT LEAST(dist, ST_distance(point, ST_ExteriorRing(pol))::numeric) INTO d;
SELECT CASE WHEN ST_Within(point,pol) THEN 1 ELSE -1 END INTO w;
SELECT ST_NumInteriorRings(pol) INTO holes;
IF holes > 0 THEN
FOR i IN 1..holes
LOOP
SELECT LEAST(d, ST_distance(point, ST_InteriorRingN(pol, i))::numeric) INTO d;
END LOOP;
END IF;
IF d < dist THEN
dist:= d;
within := w;
END IF;
END LOOP;
dist := dist * within::numeric;
RETURN dist;
END;

37
src/pg/test/expected/11_kmeans_test.out
View File

@@ -1,10 +1,43 @@
\pset format unaligned
\set ECHO all
SELECT count(DISTINCT cluster_no) as clusters from cdb_crankshaft.cdb_kmeans('select * from ppoints', 2);
-- spatial kmeans
SELECT
count(DISTINCT cluster_no) as clusters
FROM
cdb_crankshaft.cdb_kmeans('select * from ppoints', 2);
clusters
2
(1 row)
SELECT count(*) clusters from (select cdb_crankshaft.CDB_WeightedMean(the_geom, value::NUMERIC), code from ppoints group by code) p;
-- weighted mean
SELECT
count(*) clusters
FROM (
SELECT
cdb_crankshaft.CDB_WeightedMean(the_geom, value::NUMERIC),
code
FROM ppoints
GROUP BY code
) p;
clusters
52
(1 row)
-- nonspatial kmeans
SELECT
cluster_label::int in (0, 1) As cluster_label,
cluster_center::json->>'col1' As cc_col1,
cluster_center::json->>'col2' As cc_col2,
silhouettes,
inertia,
rowid
FROM cdb_crankshaft.CDB_KMeansNonspatial(
'SELECT unnest(Array[1, 1, 10, 10]) As col1, ' ||
'unnest(Array[100, 100, 2, 2]) As col2, ' ||
'unnest(Array[1, 2, 3, 4]) As cartodb_id ',
Array['col1', 'col2']::text[],
2);
cluster_label|cc_col1|cc_col2|silhouettes|inertia|rowid
t|-1.0|1.0|1.0|0.0|1
t|-1.0|1.0|1.0|0.0|2
t|1.0|-1.0|1.0|0.0|3
t|1.0|-1.0|1.0|0.0|4
(4 rows)

12
src/pg/test/expected/13_pia_test.out
View File

@@ -2,6 +2,16 @@ SET client_min_messages TO WARNING;
\set ECHO none
st_astext
-------------------------------------------
POINT(-3.67484492582767 40.4395084885993)
POINT(-3.67484492582767 40.4394914243877)
(1 row)
st_astext
------------
POINT(0 0)
(1 row)
st_astext
------------
POINT(0 0)
(1 row)

32
src/pg/test/sql/11_kmeans_test.sql
View File

@@ -1,6 +1,34 @@
\pset format unaligned
\set ECHO all
SELECT count(DISTINCT cluster_no) as clusters from cdb_crankshaft.cdb_kmeans('select * from ppoints', 2);
-- spatial kmeans
SELECT
count(DISTINCT cluster_no) as clusters
FROM
cdb_crankshaft.cdb_kmeans('select * from ppoints', 2);
SELECT count(*) clusters from (select cdb_crankshaft.CDB_WeightedMean(the_geom, value::NUMERIC), code from ppoints group by code) p;
-- weighted mean
SELECT
count(*) clusters
FROM (
SELECT
cdb_crankshaft.CDB_WeightedMean(the_geom, value::NUMERIC),
code
FROM ppoints
GROUP BY code
) p;
-- nonspatial kmeans
SELECT
cluster_label::int in (0, 1) As cluster_label,
cluster_center::json->>'col1' As cc_col1,
cluster_center::json->>'col2' As cc_col2,
silhouettes,
inertia,
rowid
FROM cdb_crankshaft.CDB_KMeansNonspatial(
'SELECT unnest(Array[1, 1, 10, 10]) As col1, ' ||
'unnest(Array[100, 100, 2, 2]) As col2, ' ||
'unnest(Array[1, 2, 3, 4]) As cartodb_id ',
Array['col1', 'col2']::text[],
2);

19
src/pg/test/sql/13_pia_test.sql
View File

@@ -5,3 +5,22 @@ with a as(
select st_geomfromtext('POLYGON((-432540.453078056 4949775.20452642,-432329.947920966 4951361.232584,-431245.028163694 4952223.31516671,-429131.071033529 4951768.00415574,-424622.07505895 4952843.13503987,-423688.327170174 4953499.20752423,-424086.294349759 4954968.38274191,-423068.388925945 4954378.63345336,-423387.653225542 4953355.67417084,-420594.869840519 4953781.00230592,-416026.095299382 4951484.06849063,-412483.018546414 4951024.5410983,-410490.399661215 4954502.24032205,-408186.197521284 4956398.91417441,-407627.262358013 4959300.94633864,-406948.770061627 4959874.85407739,-404949.583326472 4959047.74518163,-402570.908447199 4953743.46829807,-400971.358683991 4952193.11680804,-403533.488084088 4949649.89857885,-406335.177028373 4950193.19571096,-407790.456731515 4952391.46015616,-412060.672398345 4950381.2389307,-410716.93482498 4949156.7509561,-408464.162289794 4943912.8940387,-409350.599394983 4942819.84896006,-408087.791091424 4942451.6711778,-407274.045613725 4940572.4807777,-404446.196589102 4939976.71501489,-402422.964843936 4940450.3670813,-401010.654464241 4939054.8061663,-397647.247369412 4940679.80737878,-395658.413346901 4940528.84765185,-395536.852462953 4938829.79565997,-394268.923462818 4938003.7277717,-393388.720249116 4934757.80596815,-392393.301362444 4934326.71675815,-392573.527618037 4932323.40974412,-393464.640141837 4931903.10653605,-393085.597275686 4931094.7353605,-398426.261165985 4929156.87541607,-398261.174361137 4926238.00816416,-394045.059966834 4925765.18668498,-392982.960705174 4926391.81893628,-393090.272694301 4927176.84692181,-391648.240010564 4924626.06386961,-391889.914625075 4923086.14787613,-394345.177314013 4923235.086036,-395550.878718795 4917812.79243978,-399009.463978251 4912927.7157945,-398948.794855767 4911941.91010796,-398092.636652078 4911806.57392519,-401991.601817112 4911722.9204501,-406225.972607907 4914505.47286319,-411104.994569885 4912569.26941163,-412925.513522316 4913030.3608866,-414630.148884835 4914436.69169949,-414207.691417276 4919205.78028405,-418306.141109809 4917994.9580478,-424184.700779621 4918938.12432889,-426816.961458921 4923664.37379373,-420956.324227126 4923381.98014807,-420186.661267781 4924286.48693378,-420943.411166194 4926812.76394433,-419779.45457046 4928527.43466337,-419768.767899344 4930681.94459216,-421911.668097113 4930432.40620397,-423482.386112205 4933451.28047252,-427272.814773717 4934151.56473242,-427144.908678797 4939731.77191996,-428982.125554848 4940522.84445172,-428986.133056516 4942437.17281266,-431237.792396792 4947309.68284815,-432476.889648814 4947791.74800037,-432540.453078056 4949775.20452642))', 3857) as g
)
SELECT st_astext(cdb_crankshaft.CDB_PIA(g)) from a;
-- square centered on 0,0 with sides of length 2
-- expectation: point(0, 0)
WITH square AS (
SELECT 'SRID=4326;POLYGON((-1 1, 1 1, 1 -1, -1 -1, -1 1))'::geometry as g
)
SELECT ST_AsText(cdb_crankshaft.CDB_PIA(g))
FROM square;
-- MultiPolygon test
-- square centered on 0,0 with sides of length 2
-- expectation: point(0, 0)
WITH square AS (
SELECT
ST_Multi('SRID=4326;POLYGON((-1 1, 1 1, 1 -1, -1 -1, -1 1))'::geometry) as g
)
SELECT ST_AsText(cdb_crankshaft.CDB_PIA(g))
FROM square

28
src/py/crankshaft/crankshaft/analysis_data_provider.py
View File

@@ -44,8 +44,32 @@ class AnalysisDataProvider(object):
return plpy.execute(query)
@verify_data
def get_nonspatial_kmeans(self, query):
"""fetch data for non-spatial kmeans"""
def get_nonspatial_kmeans(self, params):
"""
Fetch data for non-spatial k-means.
Inputs - a dict (params) with the following keys:
colnames: a (text) list of column names (e.g.,
`['andy', 'cookie']`)
id_col: the name of the id column (e.g., `'cartodb_id'`)
subquery: the subquery for exposing the data (e.g.,
SELECT * FROM favorite_things)
Output:
A SQL query for packaging the data for consumption within
`KMeans().nonspatial`. Format will be a list of length one,
with the first element a dict with keys ('rowid', 'attr1',
'attr2', ...)
"""
agg_cols = ', '.join([
'array_agg({0}) As arr_col{1}'.format(val, idx+1)
for idx, val in enumerate(params['colnames'])
])
query = '''
SELECT {cols}, array_agg({id_col}) As rowid
FROM ({subquery}) As a
'''.format(subquery=params['subquery'],
id_col=params['id_col'],
cols=agg_cols).strip()
return plpy.execute(query)
@verify_data

81
src/py/crankshaft/crankshaft/clustering/kmeans.py
View File

@@ -30,3 +30,84 @@ class Kmeans(object):
km = KMeans(n_clusters=no_clusters, n_init=no_init)
labels = km.fit_predict(zip(xs, ys))
return zip(ids, labels)
def nonspatial(self, subquery, colnames, no_clusters=5,
standardize=True, id_col='cartodb_id'):
"""
Arguments:
query (string): A SQL query to retrieve the data required to do the
k-means clustering analysis, like so:
SELECT * FROM iris_flower_data
colnames (list): a list of the column names which contain the data
of interest, like so: ['sepal_width',
'petal_width',
'sepal_length',
'petal_length']
no_clusters (int): number of clusters (greater than zero)
id_col (string): name of the input id_column
Returns:
A list of tuples with the following columns:
cluster labels: a label for the cluster that the row belongs to
centers: center of the cluster that this row belongs to
silhouettes: silhouette measure for this value
rowid: row that these values belong to (corresponds to the value in
`id_col`)
"""
import json
from sklearn import metrics
params = {
"colnames": colnames,
"subquery": subquery,
"id_col": id_col
}
data = self.data_provider.get_nonspatial_kmeans(params)
# fill array with values for k-means clustering
if standardize:
cluster_columns = _scale_data(
_extract_columns(data))
else:
cluster_columns = _extract_columns(data)
kmeans = KMeans(n_clusters=no_clusters,
random_state=0).fit(cluster_columns)
centers = [json.dumps(dict(zip(colnames, c)))
for c in kmeans.cluster_centers_[kmeans.labels_]]
silhouettes = metrics.silhouette_samples(cluster_columns,
kmeans.labels_,
metric='sqeuclidean')
return zip(kmeans.labels_,
centers,
silhouettes,
[kmeans.inertia_] * kmeans.labels_.shape[0],
data[0]['rowid'])
# -- Preprocessing steps
def _extract_columns(data):
"""
Extract the features from the query and pack them into a NumPy array
data (list of dicts): result of the kmeans request
"""
# number of columns minus rowid column
n_cols = len(data[0]) - 1
return np.array([data[0]['arr_col{0}'.format(i+1)]
for i in xrange(n_cols)],
dtype=float).T
def _scale_data(features):
"""
Scale all input columns to center on 0 with a standard devation of 1
features (numpy matrix): features of dimension (n_features, n_samples)
"""
from sklearn.preprocessing import StandardScaler
scaler = StandardScaler()
return scaler.fit_transform(features)

2
src/py/crankshaft/requirements.txt
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

2
src/py/crankshaft/setup.py
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'],

45
src/py/crankshaft/test/test_clustering_kmeans.py
View File

@@ -2,17 +2,12 @@ import unittest
import numpy as np
# from mock_plpy import MockPlPy
# plpy = MockPlPy()
#
# import sys
# sys.modules['plpy'] = plpy
from helper import fixture_file
from crankshaft.clustering import Kmeans
from crankshaft.analysis_data_provider import AnalysisDataProvider
import crankshaft.clustering as cc
from crankshaft import random_seeds
import json
from collections import OrderedDict
@@ -24,7 +19,7 @@ class FakeDataProvider(AnalysisDataProvider):
def get_spatial_kmeans(self, query):
return self.mocked_result
def get_nonspatial_kmeans(self, query, standarize):
def get_nonspatial_kmeans(self, query):
return self.mocked_result
@@ -54,3 +49,39 @@ class KMeansTest(unittest.TestCase):
self.assertEqual(len(np.unique(labels)), 2)
self.assertEqual(len(c1), 20)
self.assertEqual(len(c2), 20)
class KMeansNonspatialTest(unittest.TestCase):
"""Testing class for k-means non-spatial"""
def setUp(self):
self.params = {"subquery": "SELECT * FROM TABLE",
"n_clusters": 5}
def test_kmeans_nonspatial(self):
"""
test for k-means non-spatial
"""
# data from:
# http://scikit-learn.org/stable/modules/generated/sklearn.cluster.KMeans.html#sklearn-cluster-kmeans
data_raw = [OrderedDict([("arr_col1", [1, 1, 1, 4, 4, 4]),
("arr_col2", [2, 4, 0, 2, 4, 0]),
("rowid", [1, 2, 3, 4, 5, 6])])]
random_seeds.set_random_seeds(1234)
kmeans = Kmeans(FakeDataProvider(data_raw))
clusters = kmeans.nonspatial('subquery', ['col1', 'col2'], 2)
cl1 = clusters[0][0]
cl2 = clusters[3][0]
for idx, val in enumerate(clusters):
if idx < 3:
self.assertEqual(val[0], cl1)
else:
self.assertEqual(val[0], cl2)
# raises exception for no data
with self.assertRaises(Exception):
kmeans = Kmeans(FakeDataProvider([]))
kmeans.nonspatial('subquery', ['col1', 'col2'], 2)