removing extra parenthesis in maxp

add option for array of ones
expose floor variable and new test file
2017-02-17 09:44:17 -05:00 · 2017-02-16 16:25:14 -05:00 · 2017-02-16 15:49:09 -05:00 · 2017-01-11 15:48:02 -05:00 · 2017-01-11 13:44:38 -05:00
18 changed files with 413 additions and 919 deletions
@@ -0,0 +1,17 @@
 -- max-p regionalization
 CREATE OR REPLACE FUNCTION
  CDB_MaxP(
      subquery TEXT,
      colnames TEXT[],
      floor_variable TEXT,
      min_size int default 1,
      initial int default 99,
      geom_col TEXT DEFAULT 'the_geom',
      id_col TEXT DEFAULT 'cartodb_id')
 RETURNS TABLE (region_class text, p_val numeric, rowid bigint)
 AS $$
  from crankshaft.clustering import MaxP
  maxp = MaxP()
  return maxp.maxp(subquery, colnames, floor_variable, floor=min_size)
 $$ LANGUAGE plpythonu;
@@ -1,26 +0,0 @@
 CREATE OR REPLACE FUNCTION
 CDB_OptimAssignments(source text,
                     drain text,
                     drain_capacity text,
                     source_production text,
                     marginal_cost text,
                     dist_matrix_query text,
                     dist_rate numeric DEFAULT 0.15,
                     dist_threshold numeric DEFAULT null)
 RETURNS table(drain_id bigint, source_id int, cost numeric, amount numeric) AS $$
 from crankshaft.optimization import Optim
 def cast_val(val):
    return float(val) if val is not None else None
 params = {'dist_rate': cast_val(dist_rate),
          'dist_threshold': cast_val(dist_threshold)}
 optim = Optim(source, drain, dist_matrix_query, drain_capacity,
              source_production, marginal_cost, **params)
 x = optim.output()
 return x
 $$ LANGUAGE plpythonu;
@@ -1,44 +0,0 @@
 -- Calculate the distance matrix using underlying road network
 -- Sample usage:
 --   select * from cdb_distancematrix('drain_table'::regclass,
 --                                    'source_table'::regclass)
 CREATE OR REPLACE FUNCTION CDB_DistanceMatrix(
    origin_table regclass,
    destination_table regclass,
    transit_mode text DEFAULT 'car'
    )
    RETURNS TABLE(origin_id bigint, destination_id bigint,
                  the_geom geometry(geometry, 4326),
                  length_km numeric, duration_sec numeric)
 AS $$
 BEGIN
    RETURN QUERY
    EXECUTE format('
        WITH pairs AS (
            SELECT
                o."cartodb_id" AS origin_id,
                d."cartodb_id" AS destination_id,
                o."the_geom" AS origin_point,
                d."the_geom" AS destination_point
            FROM
                (SELECT * FROM %I) AS o,
                (SELECT * FROM %I) AS d),
        results AS (
            SELECT
                origin_id,
                destination_id,
                (cdb_route_point_to_point(origin_point,
                                          destination_point,
                                          $1)).*
            FROM pairs)
        SELECT
            origin_id::bigint AS origin_id,
            destination_id::bigint AS destination_id,
            shape AS the_geom,
            length::numeric AS length_km,
            duration::numeric AS duration_sec
        FROM results;', origin_table, destination_table)
        USING transit_mode;
    RETURN;
 END;
 $$ LANGUAGE plpgsql;
@@ -4,4 +4,3 @@ import crankshaft.clustering
 import crankshaft.space_time_dynamics
 import crankshaft.segmentation
 import analysis_data_provider
 import crankshaft.optimization
@@ -1,11 +1,9 @@
 """class for fetching data"""
 import plpy
 import pysal_utils as pu
 import numpy as np
-class AnalysisDataProvider(object):
+
-    """Analysis providers for crankshaft functions. These rely on database
+class AnalysisDataProvider:
    access through `plpy`"""
    def get_getis(self, w_type, params):
        """fetch data for getis ord's g"""
        try:
@@ -16,7 +14,7 @@ class AnalysisDataProvider(object):
                return pu.empty_zipped_array(4)
            else:
                return result
-        except plpy.SPIError as err:
+        except plpy.SPIError, err:
            plpy.error('Analysis failed: %s' % err)
    def get_markov(self, w_type, params):
@@ -29,7 +27,7 @@ class AnalysisDataProvider(object):
                return pu.empty_zipped_array(4)
            return data
-        except plpy.SPIError as err:
+        except plpy.SPIError, err:
            plpy.error('Analysis failed: %s' % err)
    def get_moran(self, w_type, params):
@@ -42,8 +40,8 @@ class AnalysisDataProvider(object):
            if len(data) == 0:
                return pu.empty_zipped_array(2)
            return data
-        except plpy.SPIError as err:
+        except plpy.SPIError, err:
-            plpy.error('Analysis failed: %s' % err)
+            plpy.error('Analysis failed: %s' % e)
            return pu.empty_zipped_array(2)
    def get_nonspatial_kmeans(self, query):
@@ -51,7 +49,7 @@ class AnalysisDataProvider(object):
        try:
            data = plpy.execute(query)
            return data
-        except plpy.SPIError as err:
+        except plpy.SPIError, err:
            plpy.error('Analysis failed: %s' % err)
    def get_spatial_kmeans(self, params):
@@ -65,154 +63,19 @@ class AnalysisDataProvider(object):
        try:
            data = plpy.execute(query)
            return data
-        except plpy.SPIError as err:
+        except plpy.SPIError, err:
            plpy.error('Analysis failed: %s' % err)
-    def get_column(self, subquery, column, dtype=float, id_col='cartodb_id',
+    def get_maxp(self, params):
-                   condition=None):
+        """fetch data for max-p"""
        """
        Retrieve the column from the specified table from a connected
        PostgreSQL database.
        Args:
            subquery (str): subquery to retrieve column from
            column (str): column to retrieve
            dtype (type): data type in column (e.g, float, int, str)
            id_col (str, optional): Column name for index. Defaults to
                `cartodb_id`.
        Returns:
            numpy.array: column from table as a NumPy array
        """
        query = '''
            SELECT array_agg("{column}" ORDER BY "{id_col}" ASC) as col
              FROM ({subquery}) As _wrap {filter}
        '''.format(subquery=subquery,
                   column=column,
                   id_col=id_col,
                   filter='WHERE {}'.format(condition) if condition else '')
        resp = plpy.execute(query)
        return np.array(resp[0]['col'], dtype=dtype)
    def get_reduced_column(self, drain_query, capacity,
                           source_query, amount,
                           dtype=float, id_col='cartodb_id'):
        """
        Retrieve the column from the specified table from a connected
        PostgreSQL database.
        Args:
            source_query (str): source_query to retrieve column from
            column (str): column to retrieve
            dtype (type): data type in column (e.g, float, int, str)
            id_col (str, optional): Column name for index. Defaults to
                `cartodb_id`.
        Returns:
            numpy.array: column from table as a NumPy array
        """
        query = '''
            WITH cte AS (
              SELECT
                  d."{capacity}" - coalesce(s."source_claimed", 0) As
                    reduced_capacity,
                  d."{id_col}"
                FROM
                  ({drain_query}) As d
              LEFT JOIN
                  (SELECT
                      "drain_id",
                      sum("{amount}") As source_claimed
                     FROM ({source_query}) As _wrap
                   GROUP BY "drain_id") As s
                  ON
                  d."{id_col}" = s."drain_id"
              )
            SELECT
                array_agg("reduced_capacity"
                          ORDER BY "{id_col}" ASC) As col
              FROM cte
        '''.format(capacity=capacity,
                   id_col=id_col,
                   drain_query=drain_query,
                   amount=amount,
                   source_query=source_query)
        resp = plpy.execute(query)
        return np.array(resp[0]['col'], dtype=dtype)
    def get_distance_matrix(self, table, origin_ids, destination_ids):
        """Transforms a SQL table origin-destination table into a distance
        matrix.
        :param query: Table that has the data needed for building the
            distance matrix. Query should have the following columns:
            - origin_id (int)
            - destination_id (int)
            - length_km (numeric)
        :type query: str
        :param origin_ids: List of origin IDs
        :type origin_ids: list of ints
        :param destination_ids: List of origin IDs
        :type destination_ids: list of ints
        :returns: 2D array of distances from all origins to all destinations
        :rtype: numpy.array
        """
        try:
-            resp = plpy.execute('''
+            query = pu.construct_neighbor_query('queen', params)
-                SELECT "origin_id", "destination_id", "length_km"
+            data = plpy.execute(query)
                FROM (SELECT * FROM "{table}") as _wrap
            '''.format(table=table))
        except plpy.SPIError as err:
            plpy.error("Failed to build distance matrix: {}".format(err))
-        pairs = {(row['origin_id'], row['destination_id']): row['length_km']
+            if len(data) == 0:
-                 for row in resp}
+                # TODO: replace with better message in PR#157
-        distance_matrix = np.array([
+                plpy.error('No non-null valued rows')
            pairs[(origin, destination)]
            for destination in destination_ids
            for origin in origin_ids
        ])
-        return np.array(distance_matrix,
+            return data
-                        dtype=float).reshape((len(destination_ids),
+        except plpy.SPIError, err:
-                                              len(origin_ids)))
+            plpy.error('Analysis failed: %s' % err)
    def get_pairwise_distances(self, drain_query, source_query,
                               id_col='cartodb_id'):
        """Retuns the pairwise distances between row i and j for all i in
        drain_query and j in source_query
        Args:
            drain_query (str): Query that exposes the `the_geom` and
                `cartodb_id` (or what is specified in `id_col`) of the dataset
                for 'drain' locations
            source_query (str): Query that exposes the `the_geom` and
                `cartodb_id` (or what is specified in `id_col`) of the dataset
                for 'source' locations
            id_col (str, optional): Column name for table index. Defaults to
                `cartodb_id`.
        Returns:
            numpy.array: A len(s) by len(d) array of distances from source i to
                drain j
        """
        query = '''
            SELECT array_agg(ST_Distance(d."the_geom"::geography,
                                         s."the_geom"::geography) / 1000.0
                             ORDER BY d."{id_col}" ASC)  as dist
              FROM ({drain_query}) AS d, ({source_query}) AS s
            GROUP BY s."{id_col}"
            ORDER BY s."{id_col}" ASC
        '''.format(drain_query=drain_query,
                   source_query=source_query,
                   id_col=id_col)
        resp = plpy.execute(query)
        # len(s) x len(d) matrix
        return np.array([np.array(row['dist'], dtype=float)
                         for row in resp], dtype=float)
@@ -2,3 +2,4 @@
 from moran import *
 from kmeans import *
 from getis import *
 from maxp import *
@@ -38,7 +38,7 @@ class Getis:
                             ("num_ngbrs", num_ngbrs)])
        result = self.data_provider.get_getis(w_type, qvals)
-        attr_vals = pu.get_attributes(result)
+        attr_vals = pu.get_attribute(result)
        # build PySAL weight object
        weight = pu.get_weight(result, w_type, num_ngbrs)
@@ -0,0 +1,85 @@
 """
  max-p clustering
 """
 import pysal as ps
 import numpy as np
 import random
 import time
 import crankshaft.pysal_utils as pu
 from crankshaft.analysis_data_provider import AnalysisDataProvider
 class MaxP:
    def __init__(self, data_provider=None):
        if data_provider:
            self.data_provider = data_provider
        else:
            self.data_provider = AnalysisDataProvider()
    def maxp(self, subquery, colnames, floor_variable=None,
            floor=1,
            geom_col='the_geom', id_col='cartodb_id',
            ):
        """
            Inputs:
            @param subquery (text): subquery to expose the data need for the
                                    analysis. This query needs to expose all
                                    of the columns in `colnames`, `id_col`, and
                                    `geom_col`
            @param colnames (list): list of column names (as strings). This is
                                    used to calculate intra-regional homogeneity.
            @param floor_variable (text): name of column variable for the floor
            @param floor (int): the minimum bound for a variable that has to be
                                obtained in each region,
            @param geom_col (text): geometry column used for calculating the
                                    spatial neighborhood
            @param id_col (text): id column used for keeping the identity of
                                  the data
            Outputs: a list of tuples with the following columns:
            classification_id: group that the geometry belongs to
            rowid: identifier from id_col
        """
        params = {'subquery': subquery,
                  'colnames': colnames,
                  'id_col': id_col,
                  'geom_col': geom_col,
                  'floor':
                  'floor_variable':floor_variable,
                  }
        resp = self.data_provider.get_maxp(params)
        attr_vals = pu.get_attributes(resp, len(colnames))
        weight = pu.get_weight(resp, w_type='queen')
        if floor_variable == None:
            floor_variable = np.ones((weight.n, 1))
        else:
            floor_column_id = colnames.index(floor_variable)
            floor_variable = attr_vals.transpose()[floor_column_id]
        start_time = time.time()
        r = ps.Maxp(weight, attr_vals,
                    floor=floor,
                    floor_variable= floor_variable,
                    initial=10)
        # print r.regions
        cluster_classes = get_cluster_classes(weight.id_order, r.regions)
        r.inference()
        # print "elapsed time: ", time.time() - start_time
        return zip(cluster_classes, [r.pvalue] * len(weight.id_order),
                   weight.id_order)
 def get_cluster_classes(ids, clusters):
    """
    """
    cluster_classes = []
    for i in ids:
        for r_id, r in enumerate(clusters):
            if i in r:
                cluster_classes.append(r_id)
    return cluster_classes
@@ -39,7 +39,7 @@ class Moran:
        result = self.data_provider.get_moran(w_type, params)
        # collect attributes
-        attr_vals = pu.get_attributes(result)
+        attr_vals = pu.get_attribute(result, 1)
        # calculate weights
        weight = pu.get_weight(result, w_type, num_ngbrs)
@@ -68,7 +68,7 @@ class Moran:
        result = self.data_provider.get_moran(w_type, params)
-        attr_vals = pu.get_attributes(result)
+        attr_vals = pu.get_attribute(result, 1)
        weight = pu.get_weight(result, w_type, num_ngbrs)
        # calculate LISA values
@@ -96,8 +96,8 @@ class Moran:
        result = self.data_provider.get_moran(w_type, params)
        # collect attributes
-        numer = pu.get_attributes(result, 1)
+        numer = pu.get_attribute(result, 1)
-        denom = pu.get_attributes(result, 2)
+        denom = pu.get_attribute(result, 2)
        weight = pu.get_weight(result, w_type, num_ngbrs)
@@ -126,8 +126,8 @@ class Moran:
        result = self.data_provider.get_moran(w_type, params)
        # collect attributes
-        numer = pu.get_attributes(result, 1)
+        numer = pu.get_attribute(result, 1)
-        denom = pu.get_attributes(result, 2)
+        denom = pu.get_attribute(result, 2)
        weight = pu.get_weight(result, w_type, num_ngbrs)
@@ -157,8 +157,8 @@ class Moran:
        result = self.data_provider.get_moran(w_type, params)
        # collect attributes
-        attr1_vals = pu.get_attributes(result, 1)
+        attr1_vals = pu.get_attribute(result, 1)
-        attr2_vals = pu.get_attributes(result, 2)
+        attr2_vals = pu.get_attribute(result, 2)
        # create weights
        weight = pu.get_weight(result, w_type, num_ngbrs)
@@ -1 +0,0 @@
 from optim import Optim
@@ -1,301 +0,0 @@
 """optimization"""
 import plpy
 import numpy as np
 import cvxopt
 from cvxopt import solvers
 from crankshaft.analysis_data_provider import AnalysisDataProvider
 class Optim(object):
    """Linear optimization class for logistics cost minimization
    Optimization for logistics
    based on models:
      - source_amount * (marginal_cost + transport_cost * distance)
    """
    def __init__(self, source_query, drain_query, dist_matrix_table,
                 capacity_column, production_column, marginal_column,
                 **kwargs):
        # set data provider - defaults to SQL database access
        self.data_provider = kwargs.get('data_provider',
                                        AnalysisDataProvider())
        # model parameters
        self.model_params = {
            'dist_cost': kwargs.get('dist_cost', 0.15),
            'dist_threshold': kwargs.get('dist_threshold', None),
            'solver': kwargs.get('solver', 'glpk')}
        self._check_model_params()
        # database ids
        self.ids = {
            'drain_free': self.data_provider.get_column(
                drain_query,
                'cartodb_id',
                id_col='cartodb_id',
                dtype=int),
            'source_free': self.data_provider.get_column(
                source_query,
                'cartodb_id',
                dtype=int,
                condition='drain_id is null'),
            'source_fixed': self.data_provider.get_column(
                source_query,
                'cartodb_id',
                dtype=int,
                condition='drain_id is not null'),
            'drain_fixed': self.data_provider.get_column(
                source_query,
                'drain_id',
                dtype=int,
                condition='drain_id is not null'
            )}
        # model data
        self.model_data = {
            'drain_capacity': self.data_provider.get_reduced_column(
                drain_query,
                capacity_column,
                source_query,
                production_column,
                id_col='cartodb_id',
                dtype=int),
            'source_amount': self.data_provider.get_column(
                source_query,
                production_column,
                condition='drain_id is null'),
            'source_amount_fixed': self.data_provider.get_column(
                source_query,
                production_column,
                condition='drain_id is not null'),
            'marginal_cost': self.data_provider.get_column(
                drain_query,
                marginal_column),
            'distance': self.data_provider.get_distance_matrix(
                dist_matrix_table,
                self.ids['source_free'],
                self.ids['drain_free']),
            'distance_fixed': self.data_provider.get_distance_matrix(
                dist_matrix_table,
                self.ids['source_fixed'],
                self.ids['drain_fixed']
            )}
        self.model_data['cost'] = self.calc_cost()
        self.n_sources = len(self.ids['source_free'])
        self.n_drains = len(self.ids['drain_free'])
    def _check_constraints(self):
        """Check if inputs are within constraints"""
        total_capacity = self.model_data['drain_capacity'].sum()
        total_amount = self.model_data['source_amount'].sum()
        if total_amount > total_capacity:
            raise ValueError("Solution not possible. Drain capacity is "
                             "smaller than total source production.")
        elif total_capacity <= 0:
            raise ValueError("Capacity must be greater than zero")
        plpy.notice('Capacity: {total_capacity}, '
                    'Amount: {total_amount} '
                    '({perc}%)'.format(total_capacity=total_capacity,
                                       total_amount=total_amount,
                                       perc=100.0 * total_amount / total_capacity))
        return None
    def _check_model_params(self):
        """Ensure model parameters are well formed"""
        if (self.model_params['dist_threshold'] <= 0 and
                self.model_params['dist_threshold'] is not None):
            raise ValueError("`dist_threshold` must be greater than zero")
        if (self.model_params['dist_cost'] is None or
                self.model_params['dist_cost'] < 0):
            raise ValueError("`dist_cost` must be greater than zero")
        if self.model_params['solver'] not in (None, 'glpk'):
            raise ValueError("`solver` must be one of 'glpk' (default) "
                             "or None.")
        return None
    def output(self):
        """Output the calculated 'optimal' assignments if solution is not infeasible.
        :returns: List of source id/drain id pairs and the associated cost of
        transport from source to drain
        :rtype: List of tuples
        """
        # retrieve fractional assignments
        assignments = self.optim()
        # crosswalks for matrix index -> cartodb_id
        drain_id_crosswalk = {}
        for idx, cid in enumerate(self.ids['drain_free']):
            # matrix index -> cartodb_id
            drain_id_crosswalk[idx] = cid
        source_id_crosswalk = {}
        for idx, cid in enumerate(self.ids['source_free']):
            # matrix index -> cartodb_id
            source_id_crosswalk[idx] = cid
        # find non-zero entries
        source_index, drain_index = np.nonzero(assignments)
        # returns:
        #   - drain_id
        #   - source_id
        #   - cost of that pairing
        #   - amount sent via that pairing
        assigned_costs = [(
            drain_id_crosswalk[drain_index[idx]],
            source_id_crosswalk[source_val],
            self.model_data['cost'][drain_index[idx], source_val],
            round(self.model_data['source_amount'][source_val] *
                  assignments[source_val, drain_index[idx]], 6),
            False
            )
                          for idx, source_val in enumerate(source_index)]
        # Fixed vals:
        #  - self.ids['source_fixed']
        #  - self.ids['drain_fixed']
        #  -
        fixed_costs = self.fixed_values()
        # plpy.notice("FIXED COSTS: {}".format(fixed_costs))
        return assigned_costs + fixed_costs
    def fixed_values(self):
        """Return the fixed source IDs, drain IDs, costs for transport, and the
        amount that is transported.
        """
        margins = {k: val for k, val in zip(self.ids['drain_free'],
                                            self.model_data['marginal_cost'])}
        self.model_data['marginal_cost_fixed'] = [margins[d]
                                                  for d in self.ids['drain_fixed']]
        fixed_costs = self.calc_cost(source='source_amount_fixed',
                                     distance='distance_fixed',
                                     margin='marginal_cost_fixed')
        # cost = [fixed_costs[self.ids['drain_fixed'][idx], source_val]
        #         for idx, source_val in enumerate(self.ids['source_fixed'])]
        return zip(self.ids['drain_fixed'],
                   self.ids['source_fixed'],
                   [1.] * len(self.ids['source_fixed']),
                   self.model_data['source_amount_fixed'],
                   [True] * len(self.ids['drain_fixed']))
    def cost_func(self, distance, waste, marginal):
        """
        cost equation
        :param distance: distance (in km)
        :type distance: float
        :param waste: number of tons of waste. This was previously calculated
        as self.model_params['amount_per_unit'] * number of people minus the recycle_rate
        :type waste: numeric
        :param marginal: intrinsic cost per ton of a plant
        :type marginal: numeric
        :returns: cost
        :rtype: numeric
        Note: dist_cost is the cost per ton (e.g., 0.15 GBP/ton)
        """
        return waste * (marginal + self.model_params['dist_cost'] * distance)
    def calc_cost(self, source='source_amount', distance='distance',
                  margin='marginal_cost'):
        """
        Populate an d x s matrix according to the cost equation
        :returns: d x s matrix of costs from area i to plant j
        :rtype: numpy.array
        """
        costs = np.array(
            [self.cost_func(dist,
                            self.model_data[source][pair[1]],
                            self.model_data[margin][pair[0]])
             for pair, dist in np.ndenumerate(self.model_data[distance])])
        return costs.reshape(self.model_data[distance].shape)
    def optim(self):
        """solve linear optimization problem
        Equations of the form:
        minimize   c'*x     by assigning x values
        subject to G*x <= h
                   A*x = b
                   0 <= x[k] <= 1
        :returns: Fractional assignments array (of 1s and 0s) of shape c.T.
            Value at position (i, j) corresponds to the fraction of source
            `i`'s supply to drain `j`.
        :rtype: numpy.array
        """
        n_pairings = self.n_sources * self.n_drains
        # ---
        # costs
        # elements chosen to minimize sum
        cost = np.nan_to_num(self.model_data['cost'])
        cost = cvxopt.matrix(cost.ravel('F'))
        # ---
        # equality constraint variables
        # each area is serviced once
        A = cvxopt.spmatrix(1.,
                            [i // self.n_drains
                             for i in range(n_pairings)],
                            range(n_pairings), tc='d')
        b = cvxopt.matrix([1.] * self.n_sources, tc='d')
        # make nan's in cost impossible
        if np.isnan(self.model_data['distance']).any():
            i_vals, j_vals = np.where(np.isnan(self.model_data['distance']))
            for idx, i_val in enumerate(i_vals):
                i = int(i_val)
                j = int(i_val * self.n_drains + j_vals[idx])
                A[i, j] = 0
        # knock out values above distance threshold
        if self.model_params['dist_threshold']:
            j_vals, i_vals = np.where(self.model_data['distance'] >
                                      self.model_params['dist_threshold'])
            for idx, ival in enumerate(i_vals):
                A[int(ival), int(ival * self.n_drains + j_vals[idx])] = 0
        # ---
        # inequality constraint variables
        # each plant never goes over capacity
        drain_capacity = cvxopt.matrix([
            cvxopt.matrix(self.model_data['drain_capacity'], tc='d'),
            cvxopt.matrix([1.] * n_pairings, tc='d'),
            cvxopt.matrix([0.] * n_pairings, tc='d')
        ])
        # inequality maxima
        ineq_maxs = cvxopt.sparse([
            cvxopt.spmatrix(
                np.repeat(self.model_data['source_amount'], self.n_drains),
                [i % self.n_drains for i in range(n_pairings)],
                range(n_pairings), tc='d'),
            cvxopt.spmatrix(1.,
                            range(n_pairings),
                            range(n_pairings)),
            cvxopt.spmatrix(-1.,
                            range(n_pairings),
                            range(n_pairings))
            ], tc='d')
        for var in (cost, ineq_maxs, drain_capacity, A, b):
            plpy.notice('size: {}'.format(var.size))
        plpy.notice('{}, {}, {}'.format(n_pairings, self.n_sources, self.n_drains))
        # solve
        sol = solvers.lp(c=cost, G=ineq_maxs, h=drain_capacity,
                         A=A, b=b, solver=self.model_params['solver'])
        if sol['status'] != 'optimal':
            raise Exception("No solution possible: {}".format(sol))
        # NOTE: assignments needs to be shaped like self.model_data['cost'].T
        return np.array(sol['x'],
                        dtype=float)\
                   .flatten()\
                   .reshape((self.model_data['cost'].shape[1],
                             self.model_data['cost'].shape[0]))
@@ -13,10 +13,10 @@ def construct_neighbor_query(w_type, query_vals):
        @param query_vals dict: values used to construct the query
    """
-    if w_type.lower() == 'knn':
+    if w_type.lower() == 'queen':
        return knn(query_vals)
    else:
        return queen(query_vals)
    else:
        return knn(query_vals)
 # Build weight object
@@ -60,16 +60,17 @@ def query_attr_select(params):
    attr_string = ""
    template = "i.\"%(col)s\"::numeric As attr%(alias_num)s, "
-    if 'time_cols' in params:
+    if ('time_cols' in params) or ('colnames' in params):
        # if markov analysis
-        attrs = params['time_cols']
+        attrs = (params['time_cols'] if 'time_cols' in params
                 else params['colnames'])
        for idx, val in enumerate(attrs):
            attr_string += template % {"col": val, "alias_num": idx + 1}
    else:
        # if moran's analysis
        attrs = [k for k in params
-                 if k not in ('id_col', 'geom_col', 'subquery',
+                 if k not in ('id_col', 'geom_col',
                              'num_ngbrs', 'subquery')]
        for idx, val in enumerate(attrs):
@@ -100,9 +101,12 @@ def query_attr_where(params):
    attr_string = []
    template = "idx_replace.\"%s\" IS NOT NULL"
-    if 'time_cols' in params:
+    # TODO: generalize to colnames or not only?
-        # markov where clauses
+    #       this would reduce the complexity of the code here
-        attrs = params['time_cols']
+    if ('time_cols' in params) or ('colnames' in params):
        # markov and max-p where clauses
        attrs = (params['time_cols'] if 'time_cols' in params
                 else params['colnames'])
        # add values to template
        for attr in attrs:
            attr_string.append(template % attr)
@@ -111,7 +115,7 @@ def query_attr_where(params):
        # get keys
        attrs = [k for k in params
-                 if k not in ('id_col', 'geom_col', 'subquery',
+                 if k not in ('id_col', 'geom_col',
                              'num_ngbrs', 'subquery')]
        # add values to template
@@ -190,10 +194,24 @@ def queen(params):
 # to add more weight methods open a ticket or pull request
-def get_attributes(query_res, attr_num=1):
+def get_attributes(query_resp, n_cols):
    """
-        @param query_res: query results with attributes and neighbors
+        Extract the time columns and bin appropriately
    """
    return np.array([[x['attr' + str(i + 1)] for x in query_resp]
                     for i in range(n_cols)], dtype=float).transpose()
 def get_attribute(query_res, attr_num=1):
    """
        Inputs:
        @param query_res: query results with attributes and other info, of the
                          form [{'attr1': ..., 'id_col': ...},
                                {'attr1': ..., 'id_col': ...},
                                ...]
        @param attr_num: attribute number (1, 2, ...)
        Returns:
        a numpy array that represents the column in 'attr' number attr_num
    """
    return np.array([x['attr' + str(attr_num)] for x in query_res],
                    dtype=np.float)
@@ -68,7 +68,7 @@ class Markov:
        weights.transform = 'r'
        # prep time data
-        t_data = get_time_data(query_result, time_cols)
+        t_data = pu.get_attributes(query_result, len(time_cols))
        sp_markov_result = ps.Spatial_Markov(t_data,
                                             weights,
@@ -88,61 +88,15 @@ class Markov:
                                  sp_markov_result.classes[:, -1])
        # find the ups and down and overall distribution of each cell
-        trend_up, trend_down, trend, volatility = get_prob_stats(prob_dist, sp_markov_result.classes[:, -1])
+        trend_up, trend_down, trend, \
            volatility = get_prob_stats(
                prob_dist,
                sp_markov_result.classes[:, -1])
        # output the results
        return zip(trend, trend_up, trend_down, volatility, weights.id_order)
 def get_time_data(markov_data, time_cols):
    """
        Extract the time columns and bin appropriately
    """
    num_attrs = len(time_cols)
    return np.array([[x['attr' + str(i)] for x in markov_data]
                     for i in range(1, num_attrs+1)], dtype=float).transpose()
 # not currently used
 def rebin_data(time_data, num_time_per_bin):
    """
        Convert an n x l matrix into an (n/m) x l matrix where the values are
         reduced (averaged) for the intervening states:
          1 2 3 4    1.5 3.5
          5 6 7 8 -> 5.5 7.5
          9 8 7 6    8.5 6.5
          5 4 3 2    4.5 2.5
          if m = 2, the 4 x 4 matrix is transformed to a 2 x 4 matrix.
        This process effectively resamples the data at a longer time span n
         units longer than the input data.
        For cases when there is a remainder (remainder(5/3) = 2), the remaining
         two columns are binned together as the last time period, while the
         first three are binned together for the first period.
        Input:
          @param time_data n x l  ndarray: measurements of an attribute at
           different time intervals
          @param num_time_per_bin int: number of columns to average into a new
           column
        Output:
          ceil(n / m) x l ndarray of resampled time series
    """
    if time_data.shape[1] % num_time_per_bin == 0:
        # if fit is perfect, then use it
        n_max = time_data.shape[1] / num_time_per_bin
    else:
        # fit remainders into an additional column
        n_max = time_data.shape[1] / num_time_per_bin + 1
    return np.array(
      [time_data[:, num_time_per_bin * i:num_time_per_bin * (i+1)].mean(axis=1)
       for i in range(n_max)]).T
 def get_prob_dist(transition_matrix, lag_indices, unit_indices):
    """
        Given an array of transition matrices, look up the probability
@@ -0,0 +1,62 @@
 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 MaxP
 from crankshaft.analysis_data_provider import AnalysisDataProvider
 import crankshaft.clustering as cc
 from crankshaft import random_seeds
 import json
 from collections import OrderedDict
 class FakeDataProvider(AnalysisDataProvider):
    def __init__(self, fixturedata):
        self.your_maxp_data = fixturedata
    def get_maxp(self, params):
        """
          Replace this function name with the one used in your algorithm,
          and make sure to use the same function signature that is written
          for this algo in analysis_data_provider.py
        """
        return self.your_maxp_data
 class MaxPTest(unittest.TestCase):
    """Testing class for max-p regionalization"""
    def setUp(self):
        self.neighbor_data = json.loads(
          open(fixture_file('maxp.json')).read())
        self.neighbor_data = self.neighbor_data['rows']
        self.params = {"subquery": "select * from fake_table",
                       "colnames": ['population','median_hh_income'],
                       "floor_variable": 'population',
                       "floor": 260000}
    def test_maxp(self):
        """
        """
        data = [{'id': d['id'],
                 'attr1': d['attr1'],
                 'attr2':d['attr2'],
                 'neighbors': d['neighbors']} for d in self.neighbor_data]
        random_seeds.set_random_seeds(1234)
        maxp = MaxP(FakeDataProvider(data))
        regions = maxp.maxp('select * from research_nothing',['population','median_hh_income'],floor_variable='population',floor=260000)
        region_labels = [a[0] for a in regions]
        data_regionalized = zip(data,region_labels)
        for i in set(region_labels):
            sum_pop = 0
            for n in data_regionalized:
                if n[1] == i:
                    sum_pop += n[0]['attr1']
            self.assertGreaterEqual(sum_pop, 260000)
@@ -1,121 +0,0 @@
 """Unit tests for the optimizaiton module"""
 import unittest
 import numpy as np
 from crankshaft.optimization import Optim
 from crankshaft.analysis_data_provider import AnalysisDataProvider
 import cvxopt
 # suppress cvxopt GLPK messages
 cvxopt.glpk.options['msg_lev'] = 'GLP_MSG_OFF'
 class RawDataProvider(AnalysisDataProvider):
    """Raw data provider for testing purposes"""
    def __init__(self, raw_data):
        self.raw_data = raw_data
    def get_column(self, table, column, dtype=float):
        """Returns requested 'column' of data"""
        if column != 'cartodb_id':
            return np.array(self.raw_data[column], dtype=dtype)
        elif table == 'drain_table':
            return np.arange(1, len(self.raw_data['capacity_col']) + 1)
        elif table == 'source_table':
            return np.arange(1, len(self.raw_data['production_col']) + 1)
    def get_pairwise_distances(self, _source, _drain):
        """Returns pairwise distances"""
        return np.array(self.raw_data['pairwise'], dtype=float)
 class OptimTest(unittest.TestCase):
    """Testing class for Optimization module"""
    def setUp(self):
        # self.data = json.loads(
        #     open(fixture_file('optim.json')).read())
        # capacity ~ 0.01 * production given waste_per_person of 0.01
        # so capacity_col = [9, 31] / 100
        self.data = {
            'all_right': {"production_col": [10, 10, 10],
                          "capacity_col": [0.09, 0.31],
                          "marginal_col": [5, 5],
                          "pairwise": [[1, 2, 3], [3, 2, 1]]},
            'all_left': {"production_col": [10, 10, 10],
                         "capacity_col": [0.31, 0.09],
                         "marginal_col": [5, 5],
                         "pairwise": [[1, 2, 3], [3, 2, 1]]},
            '2left': {"production_col": [10, 10, 10],
                      "capacity_col": [0.21, 0.11],
                      "marginal_col": [5, 5],
                      "pairwise": [[1, 2, 3], [3, 2, 1]]},
            'infeasible': {"production_col": [10, 10, 10],
                           "capacity_col": [0.19, 0.11],
                           "marginal_col": [5, 5],
                           "pairwise": [[1, 2, 3], [3, 2, 1]]}}
        self.params = {'waste_per_person': 0.01,
                       'recycle_rate': 0.0,
                       'dist_rate': 0.15,
                       'dist_threshold': None,
                       'data_provider': None}
        self.args = ('drain_table', 'source_table',
                     'capacity_col', 'production_col',
                     'marginal_col')
        # print(self.model_data)
        # print(self.model_params)
    def test_optim_output(self):
        """Test Optim().output"""
        outputs = {'all_right': [2, 2, 2],
                   'all_left': [1, 1, 1],
                   '2left': [1, 1, 2],
                   'infeasible': None}
        for k in self.data:
            if k == 'infeasible':
                continue
            self.params['data_provider'] = RawDataProvider(self.data[k])
            optim = Optim(*self.args, **self.params)
            out_vals = optim.output()
            drain_ids = [row[0] for row in out_vals]
            print(drain_ids)
            print(k)
            self.assertTrue(drain_ids == outputs[k])
        return True
    def test_check_constraints(self):
        """Test optim._check_constraints"""
        for k in self.data:
            self.params['data_provider'] = RawDataProvider(self.data[k])
            print(k)
            try:
                optim = Optim(*self.args, **self.params)
                # pylint: disable=protected-access
                constraint_check = optim._check_constraints() is None
            except ValueError as err:
                # if infeasible, catch and say it's acceptable
                print(k)
                if k == 'infeasible':
                    constraint_check = True
                    print(constraint_check)
                else:
                    raise ValueError(err)
            self.assertTrue(constraint_check)
    # def test_check_model_params(self):
    #     """Test model param defaults are correctly formed"""
    #     for k in self.data:
    #         self.params['data_provider'] = RawDataProvider(self.data[k])
    #         optim = Optim(*self.args, **self.params)
    #         # pylint: disable=protected-access
    #         model_check = optim._check_model_params() is None
    #         self.assertTrue(model_check)
    def test_optim(self):
        """Test optim.optim method"""
        # assert False
        pass
@@ -1,8 +1,10 @@
 import unittest
 import json
 import crankshaft.pysal_utils as pu
 from crankshaft import random_seeds
 from collections import OrderedDict
 from helper import fixture_file
 class PysalUtilsTest(unittest.TestCase):
@@ -35,6 +37,8 @@ class PysalUtilsTest(unittest.TestCase):
                             "subquery": "SELECT * FROM a_list",
                             "geom_col": "the_geom",
                             "num_ngbrs": 321}
        self.neighbors_data = json.loads(
            open(fixture_file('neighbors_markov.json')).read())
    def test_query_attr_select(self):
        """Test query_attr_select"""
@@ -158,3 +162,184 @@ class PysalUtilsTest(unittest.TestCase):
        ans4 = [(None, None, None, None)]
        self.assertEqual(pu.empty_zipped_array(2), ans2)
        self.assertEqual(pu.empty_zipped_array(4), ans4)
    def test_get_attributes(self):
        """Test get_time_data"""
        import numpy as np
        data = [{'attr1': d['y1995'],
                 'attr2': d['y1996'],
                 'attr3': d['y1997'],
                 'attr4': d['y1998'],
                 'attr5': d['y1999'],
                 'attr6': d['y2000'],
                 'attr7': d['y2001'],
                 'attr8': d['y2002'],
                 'attr9': d['y2003'],
                 'attr10': d['y2004'],
                 'attr11': d['y2005'],
                 'attr12': d['y2006'],
                 'attr13': d['y2007'],
                 'attr14': d['y2008'],
                 'attr15': d['y2009']} for d in self.neighbors_data]
        result = pu.get_attributes(
            data, len(['y1995', 'y1996', 'y1997', 'y1998',
                       'y1999', 'y2000', 'y2001', 'y2002',
                       'y2003', 'y2004', 'y2005', 'y2006',
                       'y2007', 'y2008', 'y2009']))
        # expected was prepared from PySAL example:
        # f = ps.open(ps.examples.get_path("usjoin.csv"))
        # pci = np.array([f.by_col[str(y)]
        #                 for y in range(1995, 2010)]).transpose()
        # rpci = pci / (pci.mean(axis = 0))
        expected = np.array(
          [[0.87654416, 0.863147, 0.85637567, 0.84811668, 0.8446154,
            0.83271652,  0.83786314, 0.85012593, 0.85509656, 0.86416612,
            0.87119375, 0.86302631,  0.86148267, 0.86252252, 0.86746356],
           [0.9188951,  0.91757931, 0.92333258, 0.92517289, 0.92552388,
            0.90746978,  0.89830489, 0.89431991, 0.88924794, 0.89815176,
            0.91832091, 0.91706054,  0.90139505, 0.87897455, 0.86216858],
           [0.82591007, 0.82548596, 0.81989793, 0.81503235, 0.81731522,
            0.78964559,  0.80584442, 0.8084998,  0.82258551, 0.82668196,
            0.82373724, 0.81814804,  0.83675961, 0.83574199, 0.84647177],
           [1.09088176, 1.08537689, 1.08456418, 1.08415404, 1.09898841,
            1.14506948,  1.12151133, 1.11160697, 1.10888621, 1.11399806,
            1.12168029, 1.13164797,  1.12958508, 1.11371818, 1.09936775],
           [1.10731446, 1.11373944, 1.13283638, 1.14472559, 1.15910025,
            1.16898201,  1.17212488, 1.14752303, 1.11843284, 1.11024964,
            1.11943471, 1.11736468,  1.10863242, 1.09642516, 1.07762337],
           [1.42269757, 1.42118434, 1.44273502, 1.43577571, 1.44400684,
            1.44184737,  1.44782832, 1.41978227, 1.39092208, 1.4059372,
            1.40788646, 1.44052766,  1.45241216, 1.43306098, 1.4174431],
           [1.13073885, 1.13110513, 1.11074708, 1.13364636, 1.13088149,
            1.10888138,  1.11856629, 1.13062931, 1.11944984, 1.12446239,
            1.11671008, 1.10880034,  1.08401709, 1.06959206, 1.07875225],
           [1.04706124, 1.04516831, 1.04253372, 1.03239987, 1.02072545,
            0.99854316,  0.9880258,  0.99669587, 0.99327676, 1.01400905,
            1.03176742, 1.040511,  1.01749645, 0.9936394,  0.98279746],
           [0.98996986, 1.00143564, 0.99491,  1.00188408, 1.00455845,
            0.99127006,  0.97925917, 0.9683482,  0.95335147, 0.93694787,
            0.94308213, 0.92232874,  0.91284091, 0.89689833, 0.88928858],
           [0.87418391, 0.86416601, 0.84425695, 0.8404494,  0.83903044,
            0.8578708,  0.86036185, 0.86107306, 0.8500772,  0.86981998,
            0.86837929, 0.87204141,  0.86633032, 0.84946077, 0.83287146],
           [1.14196118, 1.14660262, 1.14892712, 1.14909594, 1.14436624,
            1.14450183,  1.12349752, 1.12596664, 1.12213996, 1.1119989,
            1.10257792, 1.10491258,  1.11059842, 1.10509795, 1.10020097],
           [0.97282463, 0.96700147, 0.96252588, 0.9653878,  0.96057687,
            0.95831051,  0.94480909, 0.94804195, 0.95430286, 0.94103989,
            0.92122519, 0.91010201,  0.89280392, 0.89298243, 0.89165385],
           [0.94325468, 0.96436902, 0.96455242, 0.95243009, 0.94117647,
            0.9480927,  0.93539182, 0.95388718, 0.94597005, 0.96918424,
            0.94781281, 0.93466815,  0.94281559, 0.96520315, 0.96715441],
           [0.97478408, 0.98169225, 0.98712809, 0.98474769, 0.98559897,
            0.98687073,  0.99237486, 0.98209969, 0.9877653,  0.97399471,
            0.96910087, 0.98416665,  0.98423613, 0.99823861, 0.99545704],
           [0.85570269, 0.85575915, 0.85986132, 0.85693406, 0.8538012,
            0.86191535,  0.84981451, 0.85472102, 0.84564835, 0.83998883,
            0.83478547, 0.82803648,  0.8198736,  0.82265395, 0.8399404],
           [0.87022047, 0.85996258, 0.85961813, 0.85689572, 0.83947136,
            0.82785597,  0.86008789, 0.86776298, 0.86720209, 0.8676334,
            0.89179317, 0.94202108,  0.9422231,  0.93902708, 0.94479184],
           [0.90134907, 0.90407738, 0.90403991, 0.90201769, 0.90399238,
            0.90906632,  0.92693339, 0.93695966, 0.94242697, 0.94338265,
            0.91981796, 0.91108804,  0.90543476, 0.91737138, 0.94793657],
           [1.1977611,  1.18222564, 1.18439158, 1.18267865, 1.19286723,
            1.20172869,  1.21328691, 1.22624778, 1.22397075, 1.23857042,
            1.24419893, 1.23929384,  1.23418676, 1.23626739, 1.26754398],
           [1.24919678, 1.25754773, 1.26991161, 1.28020651, 1.30625667,
            1.34790023,  1.34399863, 1.32575181, 1.30795492, 1.30544841,
            1.30303302, 1.32107766,  1.32936244, 1.33001241, 1.33288462],
           [1.06768004, 1.03799276, 1.03637303, 1.02768449, 1.03296093,
            1.05059016,  1.03405057, 1.02747623, 1.03162734, 0.9961416,
            0.97356208, 0.94241549,  0.92754547, 0.92549227, 0.92138102],
           [1.09475614, 1.11526796, 1.11654299, 1.13103948, 1.13143264,
            1.13889622,  1.12442212, 1.13367018, 1.13982256, 1.14029944,
            1.11979401, 1.10905389,  1.10577769, 1.11166825, 1.09985155],
           [0.76530058, 0.76612841, 0.76542451, 0.76722683, 0.76014284,
            0.74480073,  0.76098396, 0.76156903, 0.76651952, 0.76533288,
            0.78205934, 0.76842416,  0.77487118, 0.77768683, 0.78801192],
           [0.98391336, 0.98075816, 0.98295341, 0.97386015, 0.96913803,
            0.97370819,  0.96419154, 0.97209861, 0.97441313, 0.96356162,
            0.94745352, 0.93965462,  0.93069645, 0.94020973, 0.94358232],
           [0.83561828, 0.82298088, 0.81738502, 0.81748588, 0.80904801,
            0.80071489,  0.83358256, 0.83451613, 0.85175032, 0.85954307,
            0.86790024, 0.87170334,  0.87863799, 0.87497981, 0.87888675],
           [0.98845573, 1.02092428, 0.99665283, 0.99141823, 0.99386619,
            0.98733195,  0.99644997, 0.99669587, 1.02559097, 1.01116651,
            0.99988024, 0.97906749,  0.99323123, 1.00204939, 0.99602148],
           [1.14930913, 1.15241949, 1.14300962, 1.14265542, 1.13984683,
            1.08312397,  1.05192626, 1.04230892, 1.05577278, 1.08569751,
            1.12443486, 1.08891079,  1.08603695, 1.05997314, 1.02160943],
           [1.11368269, 1.1057147,  1.11893431, 1.13778669, 1.1432272,
            1.18257029,  1.16226243, 1.16009196, 1.14467789, 1.14820235,
            1.12386598, 1.12680236,  1.12357937, 1.1159258,  1.12570828],
           [1.30379431, 1.30752186, 1.31206366, 1.31532267, 1.30625667,
            1.31210239,  1.29989156, 1.29203193, 1.27183516, 1.26830786,
            1.2617743,  1.28656675,  1.29734097, 1.29390205, 1.29345446],
           [0.83953719, 0.82701448, 0.82006005, 0.81188876, 0.80294864,
            0.78772975,  0.82848011, 0.8259679,  0.82435705, 0.83108634,
            0.84373784, 0.83891093,  0.84349247, 0.85637272, 0.86539395],
           [1.23450087, 1.2426022,  1.23537935, 1.23581293, 1.24522626,
            1.2256767,  1.21126648, 1.19377804, 1.18355337, 1.19674434,
            1.21536573, 1.23653297,  1.27962009, 1.27968392, 1.25907738],
           [0.9769662,  0.97400719, 0.98035944, 0.97581531, 0.95543282,
            0.96480308,  0.94686376, 0.93679073, 0.92540049, 0.92988835,
            0.93442917, 0.92100464,  0.91475304, 0.90249622, 0.9021363],
           [0.84986886, 0.8986851,  0.84295997, 0.87280534, 0.85659368,
            0.88937573,  0.894401, 0.90448993, 0.95495898, 0.92698333,
            0.94745352, 0.92562488,  0.96635366, 1.02520312, 1.0394296],
           [1.01922808, 1.00258203, 1.00974428, 1.00303417, 0.99765073,
            1.00759019,  0.99192968, 0.99747298, 0.99550759, 0.97583768,
            0.9610168,  0.94779638,  0.93759089, 0.93353431, 0.94121705],
           [0.86367411, 0.85558932, 0.85544346, 0.85103025, 0.84336613,
            0.83434854,  0.85813595, 0.84667961, 0.84374558, 0.85951183,
            0.87194227, 0.89455097,  0.88283929, 0.90349491, 0.90600675],
           [1.00947534, 1.00411055, 1.00698819, 0.99513687, 0.99291086,
            1.00581626,  0.98850522, 0.99291168, 0.98983209, 0.97511924,
            0.96134615, 0.96382634,  0.95011401, 0.9434686,  0.94637765],
           [1.05712571, 1.05459419, 1.05753012, 1.04880786, 1.05103857,
            1.04800023,  1.03024941, 1.04200483, 1.0402554,  1.03296979,
            1.02191682, 1.02476275,  1.02347523, 1.02517684, 1.04359571],
           [1.07084189, 1.06669497, 1.07937623, 1.07387988, 1.0794043,
            1.0531801,  1.07452771, 1.09383478, 1.1052447,  1.10322136,
            1.09167939, 1.08772756,  1.08859544, 1.09177338, 1.1096083],
           [0.86719222, 0.86628896, 0.86675156, 0.86425632, 0.86511809,
            0.86287327,  0.85169796, 0.85411285, 0.84886336, 0.84517414,
            0.84843858, 0.84488343,  0.83374329, 0.82812044, 0.82878599],
           [0.88389211, 0.92288667, 0.90282398, 0.91229186, 0.92023286,
            0.92652175,  0.94278865, 0.93682452, 0.98655146, 0.992237,
            0.9798497,  0.93869677,  0.96947771, 1.00362626, 0.98102351],
           [0.97082064, 0.95320233, 0.94534081, 0.94215593, 0.93967,
            0.93092109,  0.92662519, 0.93412152, 0.93501274, 0.92879506,
            0.92110542, 0.91035556,  0.90430364, 0.89994694, 0.90073864],
           [0.95861858, 0.95774543, 0.98254811, 0.98919472, 0.98684824,
            0.98882205,  0.97662234, 0.95601578, 0.94905385, 0.94934888,
            0.97152609, 0.97163004,  0.9700702,  0.97158948, 0.95884908],
           [0.83980439, 0.84726737, 0.85747,  0.85467221, 0.8556751,
            0.84818516,  0.85265681, 0.84502402, 0.82645665, 0.81743586,
            0.83550406, 0.83338919,  0.83511679, 0.82136617, 0.80921874],
           [0.95118156, 0.9466212,  0.94688098, 0.9508583,  0.9512441,
            0.95440787,  0.96364363, 0.96804412, 0.97136214, 0.97583768,
            0.95571724, 0.96895368,  0.97001634, 0.97082733, 0.98782366],
           [1.08910044, 1.08248968, 1.08492895, 1.08656923, 1.09454249,
            1.10558188,  1.1214086,  1.12292577, 1.13021031, 1.13342735,
            1.14686068, 1.14502975,  1.14474747, 1.14084037, 1.16142926],
           [1.06336033, 1.07365823, 1.08691496, 1.09764846, 1.11669863,
            1.11856702,  1.09764283, 1.08815849, 1.08044313, 1.09278827,
            1.07003204, 1.08398066,  1.09831768, 1.09298232, 1.09176125],
           [0.79772065, 0.78829196, 0.78581151, 0.77615922, 0.77035744,
            0.77751194,  0.79902974, 0.81437881, 0.80788828, 0.79603865,
            0.78966436, 0.79949807,  0.80172182, 0.82168155, 0.85587911],
           [1.0052447,  1.00007696, 1.00475899, 1.00613942, 1.00639561,
            1.00162979,  0.99860739, 1.00814981, 1.00574316, 0.99030032,
            0.97682565, 0.97292596,  0.96519561, 0.96173403, 0.95890284],
           [0.95808419, 0.9382568,  0.9654441,  0.95561201, 0.96987289,
            0.96608031,  0.99727185, 1.00781194, 1.03484236, 1.05333619,
            1.0983263,  1.1704974,  1.17025154, 1.18730553, 1.14242645]])
        self.assertTrue(np.allclose(result, expected))
        self.assertTrue(type(result) == type(expected))
        self.assertTrue(result.shape == expected.shape)
@@ -105,204 +105,6 @@ class SpaceTimeTests(unittest.TestCase):
             ) in zip(result, expected):
            self.assertAlmostEqual(res_trend, exp_trend)
    def test_get_time_data(self):
        """Test get_time_data"""
        data = [{'attr1': d['y1995'],
                 'attr2': d['y1996'],
                 'attr3': d['y1997'],
                 'attr4': d['y1998'],
                 'attr5': d['y1999'],
                 'attr6': d['y2000'],
                 'attr7': d['y2001'],
                 'attr8': d['y2002'],
                 'attr9': d['y2003'],
                 'attr10': d['y2004'],
                 'attr11': d['y2005'],
                 'attr12': d['y2006'],
                 'attr13': d['y2007'],
                 'attr14': d['y2008'],
                 'attr15': d['y2009']} for d in self.neighbors_data]
        result = std.get_time_data(data, ['y1995', 'y1996', 'y1997', 'y1998',
                                          'y1999', 'y2000', 'y2001', 'y2002',
                                          'y2003', 'y2004', 'y2005', 'y2006',
                                          'y2007', 'y2008', 'y2009'])
        # expected was prepared from PySAL example:
        # f = ps.open(ps.examples.get_path("usjoin.csv"))
        # pci = np.array([f.by_col[str(y)]
        #                 for y in range(1995, 2010)]).transpose()
        # rpci = pci / (pci.mean(axis = 0))
        expected = np.array(
          [[0.87654416, 0.863147, 0.85637567, 0.84811668, 0.8446154,
            0.83271652,  0.83786314, 0.85012593, 0.85509656, 0.86416612,
            0.87119375, 0.86302631,  0.86148267, 0.86252252, 0.86746356],
           [0.9188951,  0.91757931, 0.92333258, 0.92517289, 0.92552388,
            0.90746978,  0.89830489, 0.89431991, 0.88924794, 0.89815176,
            0.91832091, 0.91706054,  0.90139505, 0.87897455, 0.86216858],
           [0.82591007, 0.82548596, 0.81989793, 0.81503235, 0.81731522,
            0.78964559,  0.80584442, 0.8084998,  0.82258551, 0.82668196,
            0.82373724, 0.81814804,  0.83675961, 0.83574199, 0.84647177],
           [1.09088176, 1.08537689, 1.08456418, 1.08415404, 1.09898841,
            1.14506948,  1.12151133, 1.11160697, 1.10888621, 1.11399806,
            1.12168029, 1.13164797,  1.12958508, 1.11371818, 1.09936775],
           [1.10731446, 1.11373944, 1.13283638, 1.14472559, 1.15910025,
            1.16898201,  1.17212488, 1.14752303, 1.11843284, 1.11024964,
            1.11943471, 1.11736468,  1.10863242, 1.09642516, 1.07762337],
           [1.42269757, 1.42118434, 1.44273502, 1.43577571, 1.44400684,
            1.44184737,  1.44782832, 1.41978227, 1.39092208, 1.4059372,
            1.40788646, 1.44052766,  1.45241216, 1.43306098, 1.4174431],
           [1.13073885, 1.13110513, 1.11074708, 1.13364636, 1.13088149,
            1.10888138,  1.11856629, 1.13062931, 1.11944984, 1.12446239,
            1.11671008, 1.10880034,  1.08401709, 1.06959206, 1.07875225],
           [1.04706124, 1.04516831, 1.04253372, 1.03239987, 1.02072545,
            0.99854316,  0.9880258,  0.99669587, 0.99327676, 1.01400905,
            1.03176742, 1.040511,  1.01749645, 0.9936394,  0.98279746],
           [0.98996986, 1.00143564, 0.99491,  1.00188408, 1.00455845,
            0.99127006,  0.97925917, 0.9683482,  0.95335147, 0.93694787,
            0.94308213, 0.92232874,  0.91284091, 0.89689833, 0.88928858],
           [0.87418391, 0.86416601, 0.84425695, 0.8404494,  0.83903044,
            0.8578708,  0.86036185, 0.86107306, 0.8500772,  0.86981998,
            0.86837929, 0.87204141,  0.86633032, 0.84946077, 0.83287146],
           [1.14196118, 1.14660262, 1.14892712, 1.14909594, 1.14436624,
            1.14450183,  1.12349752, 1.12596664, 1.12213996, 1.1119989,
            1.10257792, 1.10491258,  1.11059842, 1.10509795, 1.10020097],
           [0.97282463, 0.96700147, 0.96252588, 0.9653878,  0.96057687,
            0.95831051,  0.94480909, 0.94804195, 0.95430286, 0.94103989,
            0.92122519, 0.91010201,  0.89280392, 0.89298243, 0.89165385],
           [0.94325468, 0.96436902, 0.96455242, 0.95243009, 0.94117647,
            0.9480927,  0.93539182, 0.95388718, 0.94597005, 0.96918424,
            0.94781281, 0.93466815,  0.94281559, 0.96520315, 0.96715441],
           [0.97478408, 0.98169225, 0.98712809, 0.98474769, 0.98559897,
            0.98687073,  0.99237486, 0.98209969, 0.9877653,  0.97399471,
            0.96910087, 0.98416665,  0.98423613, 0.99823861, 0.99545704],
           [0.85570269, 0.85575915, 0.85986132, 0.85693406, 0.8538012,
            0.86191535,  0.84981451, 0.85472102, 0.84564835, 0.83998883,
            0.83478547, 0.82803648,  0.8198736,  0.82265395, 0.8399404],
           [0.87022047, 0.85996258, 0.85961813, 0.85689572, 0.83947136,
            0.82785597,  0.86008789, 0.86776298, 0.86720209, 0.8676334,
            0.89179317, 0.94202108,  0.9422231,  0.93902708, 0.94479184],
           [0.90134907, 0.90407738, 0.90403991, 0.90201769, 0.90399238,
            0.90906632,  0.92693339, 0.93695966, 0.94242697, 0.94338265,
            0.91981796, 0.91108804,  0.90543476, 0.91737138, 0.94793657],
           [1.1977611,  1.18222564, 1.18439158, 1.18267865, 1.19286723,
            1.20172869,  1.21328691, 1.22624778, 1.22397075, 1.23857042,
            1.24419893, 1.23929384,  1.23418676, 1.23626739, 1.26754398],
           [1.24919678, 1.25754773, 1.26991161, 1.28020651, 1.30625667,
            1.34790023,  1.34399863, 1.32575181, 1.30795492, 1.30544841,
            1.30303302, 1.32107766,  1.32936244, 1.33001241, 1.33288462],
           [1.06768004, 1.03799276, 1.03637303, 1.02768449, 1.03296093,
            1.05059016,  1.03405057, 1.02747623, 1.03162734, 0.9961416,
            0.97356208, 0.94241549,  0.92754547, 0.92549227, 0.92138102],
           [1.09475614, 1.11526796, 1.11654299, 1.13103948, 1.13143264,
            1.13889622,  1.12442212, 1.13367018, 1.13982256, 1.14029944,
            1.11979401, 1.10905389,  1.10577769, 1.11166825, 1.09985155],
           [0.76530058, 0.76612841, 0.76542451, 0.76722683, 0.76014284,
            0.74480073,  0.76098396, 0.76156903, 0.76651952, 0.76533288,
            0.78205934, 0.76842416,  0.77487118, 0.77768683, 0.78801192],
           [0.98391336, 0.98075816, 0.98295341, 0.97386015, 0.96913803,
            0.97370819,  0.96419154, 0.97209861, 0.97441313, 0.96356162,
            0.94745352, 0.93965462,  0.93069645, 0.94020973, 0.94358232],
           [0.83561828, 0.82298088, 0.81738502, 0.81748588, 0.80904801,
            0.80071489,  0.83358256, 0.83451613, 0.85175032, 0.85954307,
            0.86790024, 0.87170334,  0.87863799, 0.87497981, 0.87888675],
           [0.98845573, 1.02092428, 0.99665283, 0.99141823, 0.99386619,
            0.98733195,  0.99644997, 0.99669587, 1.02559097, 1.01116651,
            0.99988024, 0.97906749,  0.99323123, 1.00204939, 0.99602148],
           [1.14930913, 1.15241949, 1.14300962, 1.14265542, 1.13984683,
            1.08312397,  1.05192626, 1.04230892, 1.05577278, 1.08569751,
            1.12443486, 1.08891079,  1.08603695, 1.05997314, 1.02160943],
           [1.11368269, 1.1057147,  1.11893431, 1.13778669, 1.1432272,
            1.18257029,  1.16226243, 1.16009196, 1.14467789, 1.14820235,
            1.12386598, 1.12680236,  1.12357937, 1.1159258,  1.12570828],
           [1.30379431, 1.30752186, 1.31206366, 1.31532267, 1.30625667,
            1.31210239,  1.29989156, 1.29203193, 1.27183516, 1.26830786,
            1.2617743,  1.28656675,  1.29734097, 1.29390205, 1.29345446],
           [0.83953719, 0.82701448, 0.82006005, 0.81188876, 0.80294864,
            0.78772975,  0.82848011, 0.8259679,  0.82435705, 0.83108634,
            0.84373784, 0.83891093,  0.84349247, 0.85637272, 0.86539395],
           [1.23450087, 1.2426022,  1.23537935, 1.23581293, 1.24522626,
            1.2256767,  1.21126648, 1.19377804, 1.18355337, 1.19674434,
            1.21536573, 1.23653297,  1.27962009, 1.27968392, 1.25907738],
           [0.9769662,  0.97400719, 0.98035944, 0.97581531, 0.95543282,
            0.96480308,  0.94686376, 0.93679073, 0.92540049, 0.92988835,
            0.93442917, 0.92100464,  0.91475304, 0.90249622, 0.9021363],
           [0.84986886, 0.8986851,  0.84295997, 0.87280534, 0.85659368,
            0.88937573,  0.894401, 0.90448993, 0.95495898, 0.92698333,
            0.94745352, 0.92562488,  0.96635366, 1.02520312, 1.0394296],
           [1.01922808, 1.00258203, 1.00974428, 1.00303417, 0.99765073,
            1.00759019,  0.99192968, 0.99747298, 0.99550759, 0.97583768,
            0.9610168,  0.94779638,  0.93759089, 0.93353431, 0.94121705],
           [0.86367411, 0.85558932, 0.85544346, 0.85103025, 0.84336613,
            0.83434854,  0.85813595, 0.84667961, 0.84374558, 0.85951183,
            0.87194227, 0.89455097,  0.88283929, 0.90349491, 0.90600675],
           [1.00947534, 1.00411055, 1.00698819, 0.99513687, 0.99291086,
            1.00581626,  0.98850522, 0.99291168, 0.98983209, 0.97511924,
            0.96134615, 0.96382634,  0.95011401, 0.9434686,  0.94637765],
           [1.05712571, 1.05459419, 1.05753012, 1.04880786, 1.05103857,
            1.04800023,  1.03024941, 1.04200483, 1.0402554,  1.03296979,
            1.02191682, 1.02476275,  1.02347523, 1.02517684, 1.04359571],
           [1.07084189, 1.06669497, 1.07937623, 1.07387988, 1.0794043,
            1.0531801,  1.07452771, 1.09383478, 1.1052447,  1.10322136,
            1.09167939, 1.08772756,  1.08859544, 1.09177338, 1.1096083],
           [0.86719222, 0.86628896, 0.86675156, 0.86425632, 0.86511809,
            0.86287327,  0.85169796, 0.85411285, 0.84886336, 0.84517414,
            0.84843858, 0.84488343,  0.83374329, 0.82812044, 0.82878599],
           [0.88389211, 0.92288667, 0.90282398, 0.91229186, 0.92023286,
            0.92652175,  0.94278865, 0.93682452, 0.98655146, 0.992237,
            0.9798497,  0.93869677,  0.96947771, 1.00362626, 0.98102351],
           [0.97082064, 0.95320233, 0.94534081, 0.94215593, 0.93967,
            0.93092109,  0.92662519, 0.93412152, 0.93501274, 0.92879506,
            0.92110542, 0.91035556,  0.90430364, 0.89994694, 0.90073864],
           [0.95861858, 0.95774543, 0.98254811, 0.98919472, 0.98684824,
            0.98882205,  0.97662234, 0.95601578, 0.94905385, 0.94934888,
            0.97152609, 0.97163004,  0.9700702,  0.97158948, 0.95884908],
           [0.83980439, 0.84726737, 0.85747,  0.85467221, 0.8556751,
            0.84818516,  0.85265681, 0.84502402, 0.82645665, 0.81743586,
            0.83550406, 0.83338919,  0.83511679, 0.82136617, 0.80921874],
           [0.95118156, 0.9466212,  0.94688098, 0.9508583,  0.9512441,
            0.95440787,  0.96364363, 0.96804412, 0.97136214, 0.97583768,
            0.95571724, 0.96895368,  0.97001634, 0.97082733, 0.98782366],
           [1.08910044, 1.08248968, 1.08492895, 1.08656923, 1.09454249,
            1.10558188,  1.1214086,  1.12292577, 1.13021031, 1.13342735,
            1.14686068, 1.14502975,  1.14474747, 1.14084037, 1.16142926],
           [1.06336033, 1.07365823, 1.08691496, 1.09764846, 1.11669863,
            1.11856702,  1.09764283, 1.08815849, 1.08044313, 1.09278827,
            1.07003204, 1.08398066,  1.09831768, 1.09298232, 1.09176125],
           [0.79772065, 0.78829196, 0.78581151, 0.77615922, 0.77035744,
            0.77751194,  0.79902974, 0.81437881, 0.80788828, 0.79603865,
            0.78966436, 0.79949807,  0.80172182, 0.82168155, 0.85587911],
           [1.0052447,  1.00007696, 1.00475899, 1.00613942, 1.00639561,
            1.00162979,  0.99860739, 1.00814981, 1.00574316, 0.99030032,
            0.97682565, 0.97292596,  0.96519561, 0.96173403, 0.95890284],
           [0.95808419, 0.9382568,  0.9654441,  0.95561201, 0.96987289,
            0.96608031,  0.99727185, 1.00781194, 1.03484236, 1.05333619,
            1.0983263,  1.1704974,  1.17025154, 1.18730553, 1.14242645]])
        self.assertTrue(np.allclose(result, expected))
        self.assertTrue(type(result) == type(expected))
        self.assertTrue(result.shape == expected.shape)
    def test_rebin_data(self):
        """Test rebin_data"""
        # sample in double the time (even case since 10 % 2 = 0):
        #   (0+1)/2, (2+3)/2, (4+5)/2, (6+7)/2, (8+9)/2
        # = 0.5,     2.5,     4.5,     6.5,     8.5
        ans_even = np.array([(i + 0.5) * np.ones(10, dtype=float)
                             for i in range(0, 10, 2)]).T
        self.assertTrue(
          np.array_equal(std.rebin_data(self.time_data, 2), ans_even))
        # sample in triple the time (uneven since 10 % 3 = 1):
        #   (0+1+2)/3, (3+4+5)/3, (6+7+8)/3, (9)/1
        # = 1,         4,         7,         9
        ans_odd = np.array([i * np.ones(10, dtype=float)
                            for i in (1, 4, 7, 9)]).T
        self.assertTrue(
          np.array_equal(std.rebin_data(self.time_data, 3), ans_odd))
    def test_get_prob_dist(self):
        """Test get_prob_dist"""
        lag_indices = np.array([1, 2, 3, 4])
Author	SHA1	Message	Date
Michelle Ho	691ce95638	removing extra parenthesis in maxp	2017-02-17 09:44:17 -05:00
Michelle Ho	f8db8d9ba4	add option for array of ones	2017-02-16 16:25:14 -05:00
Michelle Ho	c6eb5d9d19	expose floor variable and new test file	2017-02-16 15:49:09 -05:00
Andy Eschbacher	ae767fc903	adds maxp algorithm	2017-01-11 15:48:02 -05:00
Andy Eschbacher	a9bef6ba1d	moves get_time_data to pysal_utils as get_attributes	2017-01-11 13:44:38 -05:00