<commit_before>#! /usr/bin/python2
# vim: set fileencoding=utf-8
"""Implement exact algorithm of the paper: Spatial Scan Statistics:
Approximations and Performance Study."""
import heapq
import numpy as np
import matplotlib.cm as cm
import matplotlib.colors as mcolor
from math import log
from more_query import SF_BBOX, KARTO_CONFIG, FIRST_TIME, LAST_TIME
BBOX = SF_BBOX
CSS = '#{} {{fill: {}; opacity: 0.5; stroke: {}; stroke-width: 0.5px;}}'
from more_query import k_split_bbox, bbox_to_polygon, compute_frequency, clock
from utils import to_css_hex
from shapely.geometry import shape, mapping, Point
from shapely import speedups
if speedups.available:
    speedups.enable()
import fiona
import json
import scipy.io as sio
from ProgressBar import AnimatedProgressBar
import persistent
from os.path import join as mkpath
import pymongo
from multiprocessing import Pool
ALLD = []


def get_discrepancy_function(total_m, total_b, support):
    """Return a binary function computing the Kulldorff discrepancy, given
    that there is total_m measured data, total_b background data and we want
    at least support points (otherwise return None)."""
    def discrepancy(m, b):
        global Reject
        """Compute d(m, b) or return None if it lacks support."""
        assert m <= total_m, "common!"
        if m < support or b < support:
            Reject += 1
            return None
        m_ratio = 1.0*m/total_m
        b_ratio = 1.0*b/total_b
        first_term = m_ratio*log(m_ratio/b_ratio)
        if abs(1-m_ratio) < 1e-7:
            # because 0 log 0 → 0 (np.allclose is probably better but 1000
            # times slower)
            second_term = 0
        else:
            second_term = (1-m_ratio)*log((1-m_ratio)/(1-b_ratio))
        if m_ratio*(1-b_ratio) > b_ratio*(1-m_ratio):
            return first_term + second_term
        else:
            return None
    return discrepancy


def add_maybe(new_value, values_so_far, max_nb_values):
    """Consider adding new_value to values_so_far if it is one of the largest
    max_nb_values."""
    real_value = new_value[0]
    if real_value is None:
        return values_so_far
    bottom_box = index_to_rect(new_value[1][0])
    upper_box = index_to_rect(new_value[1][1])
    total_box = bottom_box[:2] + upper_box[2:]
    poly = shape(bbox_to_polygon(total_box, False))
    new_value = (real_value, poly, new_value[2], new_value[3])
    ALLD.append(real_value)
    if len(values_so_far) == 0:
        return [new_value]
    if real_value > values_so_far[0][0]:
        # need_heapify = False
        # discard = False
        # for i, previous in enumerate(values_so_far):
        #     if poly.intersects(previous[1]) and not poly.touches(previous[1]):
        #         if poly.area > previous[1].area:
        #             del values_so_far[i]
        #             print('\ndelete\n')
        #             need_heapify = True
        #         else:
        #             discard = True
        #             break
        # if discard:
        #     return [new_value]
        # if not all([(previous[1].disjoint(poly) or previous[1].touches(poly))
        #             for previous in values_so_far]):
        #     return values_so_far
        # if need_heapify:
        #     heapq.heapify(values_so_far)
        if len(values_so_far) < max_nb_values:
            heapq.heappush(values_so_far, new_value)
        else:
            heapq.heappushpop(values_so_far, new_value)
    return values_so_far


def exact_grid(measured, background, discrepancy, nb_loc=1, max_size=5):
    """Given the two g×g arrays representing the measure of interest and the
    background data, find the nb_loc region that have the most discrepancy
    according to the provided binary function to compute it. Consider only<commit_msg>grow top region even when the maximum is found first<commit_after>#! /usr/bin/python2
# vim: set fileencoding=utf-8
"""Implement exact algorithm of the paper: Spatial Scan Statistics:
Approximations and Performance Study."""
import heapq
import numpy as np
import matplotlib.cm as cm
import matplotlib.colors as mcolor
from math import log
from more_query import SF_BBOX, KARTO_CONFIG, FIRST_TIME, LAST_TIME
BBOX = SF_BBOX
CSS = '#{} {{fill: {}; opacity: 0.5; stroke: {}; stroke-width: 0.5px;}}'
from more_query import k_split_bbox, bbox_to_polygon, compute_frequency, clock
from utils import to_css_hex
from shapely.geometry import shape, mapping, Point
from shapely import speedups
if speedups.available:
    speedups.enable()
import fiona
import json
import scipy.io as sio
from ProgressBar import AnimatedProgressBar
import persistent
from os.path import join as mkpath
import pymongo
from multiprocessing import Pool
ALLD = []


def get_discrepancy_function(total_m, total_b, support):
    """Return a binary function computing the Kulldorff discrepancy, given
    that there is total_m measured data, total_b background data and we want
    at least support points (otherwise return None)."""
    def discrepancy(m, b):
        global Reject
        """Compute d(m, b) or return None if it lacks support."""
        assert m <= total_m, "common!"
        if m < support or b < support:
            Reject += 1
            return None
        m_ratio = 1.0*m/total_m
        b_ratio = 1.0*b/total_b
        first_term = m_ratio*log(m_ratio/b_ratio)
        if abs(1-m_ratio) < 1e-7:
            # because 0 log 0 → 0 (np.allclose is probably better but 1000
            # times slower)
            second_term = 0
        else:
            second_term = (1-m_ratio)*log((1-m_ratio)/(1-b_ratio))
        if m_ratio*(1-b_ratio) > b_ratio*(1-m_ratio):
            return first_term + second_term
        else:
            return None
    return discrepancy


def add_maybe(new_value, values_so_far, max_nb_values):
    """Consider adding new_value to values_so_far if it is one of the largest
    max_nb_values."""
    real_value = new_value[0]
    if real_value is None:
        return values_so_far
    bottom_box = index_to_rect(new_value[1][0])
    upper_box = index_to_rect(new_value[1][1])
    total_box = bottom_box[:2] + upper_box[2:]
    poly = shape(bbox_to_polygon(total_box, False))
    new_value = (real_value, poly, new_value[2], new_value[3])
    ALLD.append(real_value)
    if len(values_so_far) == 0:
        return [new_value]
    if len(values_so_far) < max_nb_values or real_value > values_so_far[0][0]:
        # need_heapify = False
        # discard = False
        # for i, previous in enumerate(values_so_far):
        #     if poly.intersects(previous[1]) and not poly.touches(previous[1]):
        #         if poly.area > previous[1].area:
        #             del values_so_far[i]
        #             print('\ndelete\n')
        #             need_heapify = True
        #         else:
        #             discard = True
        #             break
        # if discard:
        #     return [new_value]
        # if not all([(previous[1].disjoint(poly) or previous[1].touches(poly))
        #             for previous in values_so_far]):
        #     return values_so_far
        # if need_heapify:
        #     heapq.heapify(values_so_far)
        if len(values_so_far) < max_nb_values:
            heapq.heappush(values_so_far, new_value)
        else:
            heapq.heappushpop(values_so_far, new_value)
    return values_so_far


def exact_grid(measured, background, discrepancy, nb_loc=1, max_size=5):
    """Given the two g×g arrays representing the measure of interest and the
    background data, find the nb_loc region that have the most discrepancy
    according to the provided binary function to compute it. Consider only
{'fraction_non_alphanumeric': 0.08792497069167643, 'fraction_numerical': 0.0119838478572359, 'mean_word_length': 2.7127659574468086, 'pattern_counts': {'":': 0, '<': 14, '<?xml version=': 0, '>': 9, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '2590035', 'n_tokens_mistral': 2477, 'n_tokens_neox': 2386, 'n_words': 703}