def extract_stats(neurons, config):
    '''Extract stats from neurons'''

    stats = defaultdict(dict)
    for ns, modes in config['neurite'].items():
        for n in config['neurite_type']:
            n = _NEURITE_MAP[n]
            for mode in modes:
                stat_name = _stat_name(ns, mode)
                stat = eval_stats(nm.get(ns, neurons, neurite_type=n), mode)

                if stat is None or len(stat.shape) == 0:
                    stats[n.name][stat_name] = stat
                else:
                    assert stat.shape in ((3, ), ), \
                        'Statistic must create a 1x3 result'

                    for i, suffix in enumerate('XYZ'):
                        compound_stat_name = stat_name + '_' + suffix
                        stats[n.name][compound_stat_name] = stat[i]

    for ns, modes in config['neuron'].items():
        for mode in modes:
            stat_name = _stat_name(ns, mode)
            stats[stat_name] = eval_stats(nm.get(ns, neurons), mode)

    return stats

def test_section_radial_distances_endpoint():
    ref_sec_rad_dist = nf.section_radial_distances(NEURON)

    rad_dists = fst_get('section_radial_distances', NEURON)

    nt.eq_(len(rad_dists), 84)
    nt.ok_(np.all(rad_dists == ref_sec_rad_dist))

    nrns = [nm.load_neuron(os.path.join(SWC_PATH, f)) for
            f in ('point_soma_single_neurite.swc', 'point_soma_single_neurite2.swc')]
    pop = Population(nrns)
    rad_dist_nrns = [nm.get('section_radial_distances', nrn) for nrn in nrns]
    rad_dist_pop = nm.get('section_radial_distances', pop)
    assert_items_equal(rad_dist_pop, rad_dist_nrns)

def test_segment_radial_distances_displaced_neurite():
    nrns = [nm.load_neuron(os.path.join(SWC_PATH, f)) for
            f in ('point_soma_single_neurite.swc', 'point_soma_single_neurite2.swc')]

    pop = Population(nrns)

    rad_dist_nrns = []
    for nrn in nrns:
        rad_dist_nrns.extend( nm.get('segment_radial_distances', nrn))

    rad_dist_nrns = np.array(rad_dist_nrns)

    rad_dist_pop = nm.get('segment_radial_distances', pop)

    nt.ok_(np.alltrue(rad_dist_pop == rad_dist_nrns))

def histogram(neuron, feature, ax, bins=15, normed=True, cumulative=False):
    '''
    Plot a histogram of the selected feature for the population of neurons.
    Plots x-axis versus y-axis on a scatter|histogram|binned values plot.

    Parameters :

        neurons : neuron list

        feature : str
        The feature of interest.

        bins : int
        Number of bins for the histogram.

        cumulative : bool
        Sets cumulative histogram on.

        ax : axes object
            the axes in which the plot is taking place
    '''

    feature_values = nm.get(feature, neuron)
    # generate histogram
    ax.hist(feature_values, bins=bins, cumulative=cumulative, normed=normed)

def extract_stats(neurons, config):
    '''Extract stats from neurons'''

    stats = defaultdict(dict)
    for ns, modes in config['neurite'].iteritems():
        for n in config['neurite_type']:
            n = _NEURITE_MAP[n]
            for mode in modes:
                stat_name = _stat_name(ns, mode)
                stats[n.name][stat_name] = eval_stats(nm.get(ns, neurons, neurite_type=n), mode)
                L.debug('Stat: %s, Neurite: %s, Type: %s',
                        stat_name, n, type(stats[n.name][stat_name]))

    for ns, modes in config['neuron'].iteritems():
        for mode in modes:
            stat_name = _stat_name(ns, mode)
            stats[stat_name] = eval_stats(nm.get(ns, neurons), mode)

    return stats

def load_neurite_features(filepath):
    """Unpack relevant data into megadict"""
    stuff = defaultdict(lambda: defaultdict(list))
    nrns = nm.load_neurons(filepath)
    # unpack data into arrays
    for nrn in nrns:
        for t in NEURITES_:
            for feat in FEATURES:
                stuff[feat][str(t).split(".")[1]].extend(nm.get(feat, nrn, neurite_type=t))
    return stuff

def extract_data(data_path, feature):
    '''Loads a list of neurons, extracts feature
       and transforms the fitted distribution in the correct format.
       Returns the optimal distribution, corresponding parameters,
       minimun and maximum values.
    '''
    population = nm.load_neurons(data_path)

    feature_data = [nm.get(feature, n) for n in population]
    feature_data = list(chain(*feature_data))

    return stats.optimal_distribution(feature_data)

def test_segment_radial_distances_origin():
    origin = (-100, -200, -300)
    ref_segs = nf.segment_radial_distances(NEURON)
    ref_segs_origin = nf.segment_radial_distances(NEURON, origin=origin)

    rad_dists = fst_get('segment_radial_distances', NEURON)
    rad_dists_origin = fst_get('segment_radial_distances', NEURON, origin=origin)

    nt.ok_(np.all(rad_dists == ref_segs))
    nt.ok_(np.all(rad_dists_origin == ref_segs_origin))
    nt.ok_(np.all(rad_dists_origin != ref_segs))

    nrns = [nm.load_neuron(os.path.join(SWC_PATH, f)) for
            f in ('point_soma_single_neurite.swc', 'point_soma_single_neurite2.swc')]
    pop = Population(nrns)
    rad_dist_nrns = []
    for nrn in nrns:
        rad_dist_nrns.extend(nm.get('segment_radial_distances', nrn))

    rad_dist_nrns = np.array(rad_dist_nrns)
    rad_dist_pop = nm.get('segment_radial_distances', pop)
    assert_allclose(rad_dist_nrns, rad_dist_pop)

def test_load_neuron_mixed_tree_swc():
    nrn_mix =  utils.load_neuron(os.path.join(SWC_ORD_PATH, 'sample_mixed_tree_sections.swc'))
    nt.assert_items_equal(get('number_of_sections_per_neurite', nrn_mix), [5, 3])

    nt.assert_items_equal(get('number_of_sections_per_neurite', nrn_mix),
                          get('number_of_sections_per_neurite', SWC_ORD_REF))

    nt.assert_items_equal(get('number_of_segments', nrn_mix),
                          get('number_of_segments', SWC_ORD_REF))

    nt.assert_items_equal(get('total_length', nrn_mix),
                          get('total_length', SWC_ORD_REF))

def extract_data(neurons, feature, params=None):
    '''Extracts feature from a list of neurons
       and transforms the fitted distribution in the correct format.
       Returns the optimal distribution and corresponding parameters.
       Normal distribution params (mean, std)
       Exponential distribution params (loc, scale)
       Uniform distribution params (min, range)
    '''
    if params is None:
        params = {}

    feature_data = [get(FEATURE_MAP[feature], n, **params) for n in neurons]
    feature_data = list(chain(*feature_data))
    return stats.optimal_distribution(feature_data)

def test_load_neuron_section_order_break_swc():
    nrn_mix =  utils.load_neuron(os.path.join(SWC_ORD_PATH, 'sample_disordered.swc'))

    assert_items_equal(get('number_of_sections_per_neurite', nrn_mix), [5, 3])

    assert_items_equal(get('number_of_sections_per_neurite', nrn_mix),
                       get('number_of_sections_per_neurite', SWC_ORD_REF))

    assert_items_equal(get('number_of_segments', nrn_mix),
                       get('number_of_segments', SWC_ORD_REF))

    assert_items_equal(get('total_length', nrn_mix),
                       get('total_length', SWC_ORD_REF))

def test_multiple_distr(filepath):
    '''Runs the distribution fit for multiple distributions and returns
       the optimal distribution along with the corresponding parameters.
    '''
    #  load a neuron from an SWC file
    population = nm.load_neurons(filepath)

    # Create a list of basic distributions
    distr_to_check = ('norm', 'expon', 'uniform')

    # Get the soma radii of a population of neurons
    soma_size = nm.get('soma_radii', population)

    # Find the best fit distribution
    return st.optimal_distribution(soma_size, distr_to_check)

 def test_neuron_not_corrupted(self):
     # Regression for #492: dendrogram was corrupting
     # neuron used to construct it.
     # This caused the section path distance calculation
     # to raise a KeyError exception.
     get('section_path_distances', NEURON)

 def time_number_of_sections_per_neurite(self):
     nm.get('number_of_sections_per_neurite', self.neuron)

 def time_segment_taper_rates(self):
     nm.get('segment_taper_rates', self.neuron)

 def time_segment_midpoints(self):
     nm.get('segment_midpoints', self.neuron)

 def time_segment_radii(self):
     nm.get('segment_radii', self.neuron)

 def time_number_of_forking_points(self):
     nm.get('number_of_forking_points', self.neuron)

 def time_partition(self):
     nm.get('partition', self.neuron)

 def time_local_bifurcation_angles(self):
     nm.get('local_bifurcation_angles', self.neuron)