3
\$\begingroup\$

When performing a chi-squared test, one takes the square of the differences of the expected counts per bin and observed counts per bin, and divides these per-bin differences by the expected counts per bin, as seen in the formula below.

enter image description here

However, I've heard that it is good practice to modify the bins such that the observed counts are above a threshold (typically 5, sometimes less) as a large number of bin counts below such a threshold can result in a bad fit (assuming minimized chi-squared). If the observed bin count is less than this threshold, then the bin is merged with the next bin. Assuming a distribution (such as a Gaussian) with a central peak, the next bin would be the next-right bin (i to i+1 bin) when the bins are left of the central peak while the next bin would be the next-left bin (i to i-1 bin) when the bins are to the right of the central peak.

I've created an algorithm that I believe works and covers all-edge cases (assuming a single central peak). I was wondering how it could be improved in terms of speed/efficiency. I also feel like I am duplicating code using similar approaches in two while-loops; can this be averted?

import numpy as np

## CONVENIENCE FUNCTIONS
get_midpoints = lambda edges : (edges[1:] + edges[:-1])/2
verify = lambda string, data : print("\n .. {} {}:\n{}\n".format(len(data), string, data))

## SAMPLE DATA
bin_edges = np.linspace(0, 40, 9)
bin_mids = get_midpoints(bin_edges)
bin_counts = np.array([2, 12, 21, 31, 18, 9, 3, 1])

## DESIRED OUTPUT (VIA THRESHOLD=5)
upd_edges = np.array([0, 10, 15, 20, 25, 40])
upd_mids = get_midpoints(upd_edges)
upd_counts = np.array([14, 21, 31, 18, 13])

## ALGORITHM TO MODIFY BINS
def modify_bins(edges, counts, threshold=None):
    """
    This function returns the bin edges, bin midpoints, and bin counts
    with the option of including a count threshold.

    edges       :   type <array>
    counts      :   type <array>
    threshold   :   type <int> greater than 0 or None
    """
    if threshold is None:
        return edges, mids, counts
    else:
        ## INITIALIZE
        loc_peak = np.argmax(counts)
        left_edges = [edges[0]]
        left_counts = []
        right_edges = []
        right_counts = []

        ## CONSOLIDATE BINS LEFT OF PEAK
        idx = 0
        while idx <= loc_peak:
            count_tmp = counts[idx]
            edge_tmp = edges[idx+1]
            while count_tmp < threshold:
                idx += 1
                count_tmp += counts[idx]
                edge_tmp = edges[idx+1]
            left_edges.append(edge_tmp)
            left_counts.append(count_tmp)
            idx += 1

        ## CONSOLIDATE BINS RIGHT OF PEAK
        idx = len(counts) - 1
        while idx > loc_peak:
            count_tmp = counts[idx]
            edge_tmp = edges[idx+1]
            while count_tmp < threshold:
                idx -= 1
                count_tmp += counts[idx]
                edge_tmp = edges[idx+1]
            right_edges.append(edge_tmp)
            right_counts.append(count_tmp)
            idx -= 1

        ## GROUP LEFT & RIGHT DATA
        right_edges[0] = edges[-1]
        mod_edges = np.array(left_edges + right_edges[::-1])
        mod_counts = np.array(left_counts + right_counts[::-1])
        mod_mids = get_midpoints(mod_edges)
        return mod_edges, mod_mids, mod_counts

# mod_edges, mod_mids, mod_counts = modify_bins(bin_edges, bin_counts, threshold=None)
mod_edges, mod_mids, mod_counts = modify_bins(bin_edges, bin_counts, threshold=5)

## VERIFY
verify("ORIGINAL EDGES", bin_edges)
verify("ORIGINAL MIDPOINTS", bin_mids)
verify("ORIGINAL COUNTS", bin_counts)
print("-----------------------------")
verify("DESIRED EDGES", upd_edges)
verify("DESIRED MIDPOINTS", upd_mids)
verify("DESIRED COUNTS", upd_counts)
print("-----------------------------")
verify("MODIFIED EDGES", mod_edges)
verify("MODIFIED MIDPOINTS", mod_mids)
verify("MODIFIED COUNTS", mod_counts)
\$\endgroup\$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.