# -*- coding: utf-8 -*-
"""neighbourhoods.py

Scripts to extract the areal units where the different classes are
over-represented, and cluster the areal units that have common boundaries.
"""
import math

import marble as mb
from common import (regroup_per_class,
                   return_categories,
                   compute_totals)


__all__ = ["overrepresented_units",
           "neighbourhoods"]


#
# Helper functions
#


#
# Callable functions
#
def overrepresented_units(distribution, classes=None):
    """ Find the areal units in which each class is over-represented
   
    We say that a class `\alpha` is overrepresented in that tract `t` if the
    representation `r_\alpha(t)` is such that

    .. math::
        r_\alpha(t) > 1 + 2.57 \sigma_\alpha(t)

    Parameters
    ----------

    distribution: nested dictionaries
        Number of people per class, per areal unit as given in the raw data
        (ungrouped). The dictionary must have the following formatting:
        > {areal_id: {class_id: number}}

    classes: dictionary of lists
        When the original categories need to be aggregated into different
        classes. 
        > {class: [categories belonging to this class]}
        This can be arbitrarily imposed, or computed with uncover_classes
        function of this package.

    Returns
    -------

    areal_units: dictionary of lists
        Dictionnary of classes, with the list of areal units where this class is
        overrepresented with 99% confidence.
        > {class:[list of areal units]}
    """
    # Regroup into classes if specified. Otherwise return categories indicated
    # in the data
    if not classes:
       classes = return_categories(distribution) 


    ## Compute the representation of the different classes in all areal units
    rep = mb.representation(distribution, classes)

    ## Find the tracts where classes are overrepresented
    areal_units = {cl:[au for au in rep
                          if rep[au][cl][0] > 1 + 2.57*math.sqrt(rep[au][cl][1])] 
                    for cl in classes}

    return areal_units




def neighbourhoods(distribution, areal_units, classes=None):
    """ Return the neighbourhoods where different classes gather

    Parameter
    ---------

    distribution: nested dictionaries
        Number of people per class, per areal unit as given in the raw data
        (ungrouped). The dictionary must have the following formatting:
        > {areal_id: {class_id: number}}

    areal_units: dictionnary
        Dictionnary of areal unit ids with shapely polygon object representing
        the unit's geometry as values.

    classes: dictionary of lists
        When the original categories need to be aggregated into different
        classes. 
        > {class: [categories belonging to this class]}
        This can be arbitrarily imposed, or computed with uncover_classes
        function of this package.

    Returns
    -------

    neighbourhoods: dictionary
        Dictionary of classes names with list of neighbourhoods (that are
        each represented by a list of areal unit)
        > {'class': [ [areal units in cluster i], ...]}
    """
    neighbourhoods = {}

    ## Find the areal units where classes are overrepresented
    or_units = overrepresented_units(distribution, classes)

    return neighbourhoods

@@ -196,16 +196,31 @@
 ort math
+%0Aimport shapely
 %0A%0Aimport
@@ -410,16 +410,41 @@
 urhoods%22
+,%0A           %22clustering%22
 %5D%0A%0A%0A#%0A# 
@@ -462,16 +462,719 @@
 tions%0A#%0A
+def _adjacency(areal_units)%0A    %22%22%22 Compute the adjacency matrix of areal units%0A%0A    Two areal units are said adjacent if their repective boundaries touch one%0A    another (see shapely's %60touch%60 function for more information).%0A%0A    Parameter%0A    ---------%0A%0A    areal_units: dictionary%0A        Dictionnary of areal unit ids with shapely polygon object representing%0A        the unit's geometry as values.%0A%0A    Returns%0A    -------%0A%0A    adjacency: dictionary%0A    %22%22%22%0A%0A    ## Compute adjacency list%0A    adjacency = %7Ba:%5B%5D for a in areal_units%7D%0A    for a0,a1 in itertools.permutations(areal_units, 2):%0A        if blocks%5Ba1%5D.touches(areal_units%5Ba0%5D):%0A            adjacency%5Ba0%5D.append(a1)%0A%0A    return adjacency%0A%0A%0A
 %0A%0A#%0A# Ca
@@ -2077,38 +2077,32 @@
    -------%0A%0A    
-areal_
 units: dictionar
@@ -2829,17 +2829,16 @@
 _units%0A%0A
-%0A
 %0A%0Adef ne
@@ -4071,8 +4071,1012 @@
 urhoods%0A
+%0A%0A%0Adef clustering(distribution, areal_units, classes=None):%0A    %22%22%22 Return the clustering coefficient for the different classes%0A    %0A    %5BAdd the definition of clustering here%5D%0A%0A    Parameter%0A    ---------%0A%0A    distribution: nested dictionaries%0A        Number of people per class, per areal unit as given in the raw data%0A        (ungrouped). The dictionary must have the following formatting:%0A        %3E %7Bareal_id: %7Bclass_id: number%7D%7D%0A%0A    areal_units: dictionnary%0A        Dictionnary of areal unit ids with shapely polygon object representing%0A        the unit's geometry as values.%0A%0A    classes: dictionary of lists%0A        When the original categories need to be aggregated into different%0A        classes. %0A        %3E %7Bclass: %5Bcategories belonging to this class%5D%7D%0A        This can be arbitrarily imposed, or computed with uncover_classes%0A        function of this package.%0A%0A    Returns%0A    -------%0A%0A    clustering: dictionary%0A        Dictionary of classes names with clustering values.%0A    %22%22%22%0A    return%0A

Add function to compute the adjacency matrix