# Path: modality/critical_bandwidth.py
# def critical_bandwidth_m_modes(data, m, I=(-np.inf, np.inf), htol=1e-3):
#         # I is interval over which density is tested for unimodality
#     hmax = (np.max(data)-np.min(data))/2.0
#     return bisection_search_most_m_modes(0, hmax, htol, data, m, I)


import numpy as np
    import matplotlib.pyplot as plt
from scipy.signal import argrelextrema
from sklearn.neighbors import KernelDensity
from .critical_bandwidth import critical_bandwidth_m_modes
from __future__ import unicode_literals



def best_split(data, I=(-np.inf, np.inf)):
    '''With bimodal data, finding split at lowest density.'''
    h_crit = critical_bandwidth_m_modes(data, 2, I)
    kde = KernelDensity(kernel='gaussian', bandwidth=h_crit).fit(data.reshape(-1, 1))
    x = np.linspace(max(np.min(data), I[0]), min(np.max(data), I[1]), 200)
    y = np.exp(kde.score_samples(x.reshape(-1, 1)))
    modes = argrelextrema(np.hstack([[0], y, [0]]), np.greater)[0]
    if len(modes) != 2:
        raise ValueError("{} modes at: {}".format(len(modes), x[modes-1]))
    ind_min = modes[0]-1 + argrelextrema(y[(modes[0]-1):(modes[1]-1)], np.less)[0]
    return x[ind_min]

if __name__ == '__main__':
    if 1:
        N = 1000
        data = np.hstack([np.random.randn(N/2), np.random.randn(N/4)+4])

# Path: modality/critical_bandwidth.py
# def critical_bandwidth_m_modes(data, m, I=(-np.inf, np.inf), htol=1e-3):
#         # I is interval over which density is tested for unimodality
#     hmax = (np.max(data)-np.min(data))/2.0
#     return bisection_search_most_m_modes(0, hmax, htol, data, m, I)


# Path: modality/best_split.py
import numpy as np
    import matplotlib.pyplot as plt
from scipy.signal import argrelextrema
from sklearn.neighbors import KernelDensity
from .critical_bandwidth import critical_bandwidth_m_modes

from __future__ import unicode_literals



def best_split(data, I=(-np.inf, np.inf)):
    '''With bimodal data, finding split at lowest density.'''
    h_crit = critical_bandwidth_m_modes(data, 2, I)
    kde = KernelDensity(kernel='gaussian', bandwidth=h_crit).fit(data.reshape(-1, 1))
    x = np.linspace(max(np.min(data), I[0]), min(np.max(data), I[1]), 200)
    y = np.exp(kde.score_samples(x.reshape(-1, 1)))
    modes = argrelextrema(np.hstack([[0], y, [0]]), np.greater)[0]
    if len(modes) != 2:
        raise ValueError("{} modes at: {}".format(len(modes), x[modes-1]))
    ind_min = modes[0]-1 + argrelextrema(y[(modes[0]-1):(modes[1]-1)], np.less)[0]
    return x[ind_min]

if __name__ == '__main__':
    if 1:
        N = 1000
        data = np.hstack([np.random.randn(N/2), np.random.randn(N/4)+4])

        h_crit = critical_bandwidth_m_modes(data, 2)