Python numpy.histogram2d() Examples

The following are 30 code examples of numpy.histogram2d(). You can vote up the ones you like or vote down the ones you don't like, and go to the original project or source file by following the links above each example. You may also want to check out all available functions/classes of the module numpy , or try the search function .
Example #1
Source File: discrete_model.py    From Splunking-Crime with GNU Affero General Public License v3.0 6 votes vote down vote up
def pred_table(self, threshold=.5):
        """
        Prediction table

        Parameters
        ----------
        threshold : scalar
            Number between 0 and 1. Threshold above which a prediction is
            considered 1 and below which a prediction is considered 0.

        Notes
        ------
        pred_table[i,j] refers to the number of times "i" was observed and
        the model predicted "j". Correct predictions are along the diagonal.
        """
        model = self.model
        actual = model.endog
        pred = np.array(self.predict() > threshold, dtype=float)
        bins = np.array([0, 0.5, 1])
        return np.histogram2d(actual, pred, bins=bins)[0] 
Example #2
Source File: mutual_information.py    From NeuroKit with MIT License 6 votes vote down vote up
def _mutual_information_varoquaux(x, y, bins=256, sigma=1, normalized=True):
    """Based on Gael Varoquaux's implementation: https://gist.github.com/GaelVaroquaux/ead9898bd3c973c40429."""
    jh = np.histogram2d(x, y, bins=bins)[0]

    # smooth the jh with a gaussian filter of given sigma
    scipy.ndimage.gaussian_filter(jh, sigma=sigma, mode="constant", output=jh)

    # compute marginal histograms
    jh = jh + np.finfo(float).eps
    sh = np.sum(jh)
    jh = jh / sh
    s1 = np.sum(jh, axis=0).reshape((-1, jh.shape[0]))
    s2 = np.sum(jh, axis=1).reshape((jh.shape[1], -1))

    if normalized:
        mi = ((np.sum(s1 * np.log(s1)) + np.sum(s2 * np.log(s2))) / np.sum(jh * np.log(jh))) - 1
    else:
        mi = np.sum(jh * np.log(jh)) - np.sum(s1 * np.log(s1)) - np.sum(s2 * np.log(s2))

    return mi 
Example #3
Source File: SS_det_only.py    From EXOSIMS with BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
def hist(self, nplan, xedges, yedges):
        """Returns completeness histogram for Monte Carlo simulation
        
        This function uses the inherited Planet Population module.
        
        Args:
            nplan (float):
                number of planets used
            xedges (float ndarray):
                x edge of 2d histogram (separation)
            yedges (float ndarray):
                y edge of 2d histogram (dMag)
        
        Returns:
            h (ndarray):
                2D numpy ndarray containing completeness histogram
        
        """
        
        s, dMag = self.genplans(nplan)
        # get histogram
        h, yedges, xedges = np.histogram2d(dMag, s.to('AU').value, bins=1000,
                range=[[yedges.min(), yedges.max()], [xedges.min(), xedges.max()]])
        
        return h, xedges, yedges 
Example #4
Source File: ext_dmdt.py    From feets with MIT License 6 votes vote down vote up
def fit(self, magnitude, time, dt_bins, dm_bins):
        def delta_calc(idx):
            t0 = time[idx]
            m0 = magnitude[idx]
            deltat = time[idx + 1 :] - t0
            deltam = magnitude[idx + 1 :] - m0

            deltat[np.where(deltat < 0)] *= -1
            deltam[np.where(deltat < 0)] *= -1

            return np.column_stack((deltat, deltam))

        lc_len = len(time)
        n_vals = int(0.5 * lc_len * (lc_len - 1))

        deltas = np.vstack(tuple(delta_calc(idx) for idx in range(lc_len - 1)))

        deltat = deltas[:, 0]
        deltam = deltas[:, 1]

        bins = [dt_bins, dm_bins]
        counts = np.histogram2d(deltat, deltam, bins=bins, normed=False)[0]
        result = np.fix(255.0 * counts / n_vals + 0.999).astype(int)

        return {"DMDT": result} 
Example #5
Source File: geographical.py    From typhon with MIT License 6 votes vote down vote up
def gridded_mean(lat, lon, data, grid):
    """Grid data along latitudes and longitudes

    Args:
        lat: Grid points along latitudes as 1xN dimensional array.
        lon: Grid points along longitudes as 1xM dimensional array.
        data: The data as NxM numpy array.
        grid: A tuple with two numpy arrays, consisting of latitude and
            longitude grid points.

    Returns:
        Two matrices in grid form: the mean and the number of points of `data`.
    """
    grid_sum, _, _ = np.histogram2d(lat, lon, grid, weights=data)
    grid_number, _, _ = np.histogram2d(lat, lon, grid)

    return grid_sum / grid_number, grid_number 
Example #6
Source File: OptimalProjection.py    From scattertext with Apache License 2.0 6 votes vote down vote up
def morista_index(points):
    # Morisita Index of Dispersion

    N = points.shape[1]

    ims = []
    for i in range(1, N):
        bins, _, _ = np.histogram2d(points[0], points[1], i)

        # I_M  = Q * (\sum_{k=1}^{Q}{n_k * (n_k - 1)})/(N * (N _ 1))
        Q = len(bins)  # num_quadrants

        # Eqn 1.
        I_M = Q * np.sum(np.ravel(bins) * (np.ravel(bins) - 1)) / (N * (N - 1))
        ims.append([i, I_M])


    return np.array(ims).T[1].max() 
Example #7
Source File: discrete_model.py    From vnpy_crypto with MIT License 6 votes vote down vote up
def pred_table(self, threshold=.5):
        """
        Prediction table

        Parameters
        ----------
        threshold : scalar
            Number between 0 and 1. Threshold above which a prediction is
            considered 1 and below which a prediction is considered 0.

        Notes
        ------
        pred_table[i,j] refers to the number of times "i" was observed and
        the model predicted "j". Correct predictions are along the diagonal.
        """
        model = self.model
        actual = model.endog
        pred = np.array(self.predict() > threshold, dtype=float)
        bins = np.array([0, 0.5, 1])
        return np.histogram2d(actual, pred, bins=bins)[0] 
Example #8
Source File: abstract_net.py    From bonnet with GNU General Public License v3.0 6 votes vote down vote up
def pix_histogram(self, mask, lbl):
    '''
      get individual mask and label and create 2d hist
    '''
    # flatten mask and cast
    flat_mask = mask.flatten().astype(np.uint32)
    # flatten label and cast
    flat_label = lbl.flatten().astype(np.uint32)
    # get the histogram
    histrange = np.array([[-0.5, self.num_classes - 0.5],
                          [-0.5, self.num_classes - 0.5]], dtype='float64')
    h_now, _, _ = np.histogram2d(np.array(flat_mask),
                                 np.array(flat_label),
                                 bins=self.num_classes,
                                 range=histrange)
    return h_now 
Example #9
Source File: test_samples.py    From thewalrus with Apache License 2.0 6 votes vote down vote up
def test_displaced_two_mode_state_hafnian(self, sample_func):
        """Test the sampling routines by comparing the photon number frequencies and the exact
        probability distribution of a two mode coherent state
        """
        n_samples = 1000
        n_cut = 6
        sigma = np.identity(4)
        mean = 5 * np.array([0.1, 0.25, 0.1, 0.25])
        samples = sample_func(sigma, samples=n_samples, mean=mean, cutoff=n_cut)
        # samples = hafnian_sample_classical_state(sigma, mean = mean, samples = n_samples)
        probs = np.real_if_close(
            np.array(
                [
                    [density_matrix_element(mean, sigma, [i, j], [i, j]) for i in range(n_cut)]
                    for j in range(n_cut)
                ]
            )
        )
        freq, _, _ = np.histogram2d(samples[:, 1], samples[:, 0], bins=np.arange(0, n_cut + 1))
        rel_freq = freq / n_samples

        assert np.allclose(
            rel_freq, probs, rtol=rel_tol / np.sqrt(n_samples), atol=rel_tol / np.sqrt(n_samples)
        ) 
Example #10
Source File: pid_analysis.py    From flight_review with BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
def hist2d(self, x, y, weights, bins):   #bins[nx,ny]
        ### generates a 2d hist from input 1d axis for x,y. repeats them to match shape of weights X*Y (data points)
        ### x will be 0-100%
        freqs = np.repeat(np.array([y], dtype=np.float64), len(x), axis=0)
        throts = np.repeat(np.array([x], dtype=np.float64), len(y), axis=0).transpose()
        throt_hist_avr, throt_scale_avr = np.histogram(x, 101, [0, 100])

        hist2d = np.histogram2d(throts.flatten(), freqs.flatten(),
                                range=[[0, 100], [y[0], y[-1]]],
                                bins=bins, weights=weights.flatten(), normed=False)[0].transpose()

        hist2d = np.array(abs(hist2d), dtype=np.float64)
        hist2d_norm = np.copy(hist2d)
        hist2d_norm /=  (throt_hist_avr + 1e-9)

        return {'hist2d_norm':hist2d_norm, 'hist2d':hist2d, 'throt_hist':throt_hist_avr,'throt_scale':throt_scale_avr} 
Example #11
Source File: test_gate.py    From FlowCal with MIT License 6 votes vote down vote up
def test_gate_fraction_2_error_large_gate_fraction(self):
        bins = [-0.5, 0.5, 1.5, 2.5, 3.5, 4.5]
        self.assertRaises(
            ValueError,
            FlowCal.gate.density2d,
            self.pyramid,
            bins=bins,
            gate_fraction=1.1,
            sigma=0.0,
            )

    # Test implicit gating (when values exist outside specified bins)
    #
    # The expected behavior is that density2d() should mimic
    # np.histogram2d()'s behavior: values outside the specified bins are
    # ignored (in the context of a gate function, this means they are
    # implicitly gated out). 
Example #12
Source File: test_gate.py    From FlowCal with MIT License 6 votes vote down vote up
def test_sub_bin_1_mask(self):
        bins = [0.5, 2.5, 4.5]
        np.testing.assert_array_equal(
            FlowCal.gate.density2d(
                self.slope, bins=bins, gate_fraction=13.0/30, sigma=0.0,
                full_output=True).mask,
            np.array([0,0,0,0,0,0,0,0,0,0,1,1,0,0,1,1,0,0,0,0,1,1,0,1,1,
                1,1,1,1,1], dtype=bool)
            )

    # Test bins edge case (when bin edges = values)
    #
    # Again, the expected behavior is that density2d() should mimic
    # np.histogram2d()'s behavior which will group the right-most two values
    # together in the same bin (since the last bins interval is fully closed,
    # as opposed to all other bins intervals which are half-open). 
Example #13
Source File: distributional_nbd.py    From netrd with MIT License 6 votes vote down vote up
def spectral_distribution(points, cumulative=True):
    """ 
    Returns the distribution of complex values (in r,theta-space).
    """

    points = np.array([(np.abs(z), np.angle(z)) for z in points])
    r, theta = np.split(points, 2, axis=1)

    r = np.array([logr(x, r.max()) for x in r])

    Z, R, THETA = np.histogram2d(
        x=r[:, 0],
        y=theta[:, 0],
        bins=(np.linspace(0, 1, 101), np.linspace(0, np.pi, 101)),
    )

    if cumulative:
        Z = Z.cumsum(axis=0).cumsum(axis=1)
        Z = Z / Z.max()

    return Z.flatten() 
Example #14
Source File: stats.py    From scprep with GNU General Public License v3.0 6 votes vote down vote up
def mutual_information(x, y, bins=8):
    """Mutual information score with set number of bins

    Helper function for `sklearn.metrics.mutual_info_score` that builds a
    contingency table over a set number of bins.
    Credit: `Warran Weckesser <https://stackoverflow.com/a/20505476/3996580>`_.


    Parameters
    ----------
    x : array-like, shape=[n_samples]
        Input data (feature 1)
    y : array-like, shape=[n_samples]
        Input data (feature 2)
    bins : int or array-like, (default: 8)
        Passed to np.histogram2d to calculate a contingency table.

    Returns
    -------
    mi : float
        Mutual information between x and y.

    Examples
    --------
    >>> import scprep
    >>> data = scprep.io.load_csv("my_data.csv")
    >>> mi = scprep.stats.mutual_information(data['GENE1'], data['GENE2'])
    """
    x, y = _vector_coerce_two_dense(x, y)
    c_xy = np.histogram2d(x, y, bins)[0]
    mi = metrics.mutual_info_score(None, None, contingency=c_xy)
    return mi 
Example #15
Source File: BrownCompleteness.py    From EXOSIMS with BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
def hist(self, nplan, xedges, yedges):
        """Returns completeness histogram for Monte Carlo simulation
        
        This function uses the inherited Planet Population module.
        
        Args:
            nplan (float):
                number of planets used
            xedges (float ndarray):
                x edge of 2d histogram (separation)
            yedges (float ndarray):
                y edge of 2d histogram (dMag)
        
        Returns:
            float ndarray:
                2D numpy ndarray containing completeness frequencies
        
        """
        
        s, dMag = self.genplans(nplan)
        # get histogram
        h, yedges, xedges = np.histogram2d(dMag, s.to('AU').value, bins=1000,
                range=[[yedges.min(), yedges.max()], [xedges.min(), xedges.max()]])
        
        return h, xedges, yedges 
Example #16
Source File: sct_maths.py    From spinalcordtoolbox with MIT License 6 votes vote down vote up
def mutual_information(x, y, nbins=32, normalized=False):
    """
    Compute mutual information
    :param x: 1D numpy.array : flatten data from an image
    :param y: 1D numpy.array : flatten data from an image
    :param nbins: number of bins to compute the contingency matrix (only used if normalized=False)
    :return: float non negative value : mutual information
    """
    import sklearn.metrics
    if normalized:
        mi = sklearn.metrics.normalized_mutual_info_score(x, y)
    else:
        c_xy = np.histogram2d(x, y, nbins)[0]
        mi = sklearn.metrics.mutual_info_score(None, None, contingency=c_xy)
    # mi = adjusted_mutual_info_score(None, None, contingency=c_xy)
    return mi 
Example #17
Source File: mutual_information_matrix.py    From netrd with MIT License 5 votes vote down vote up
def find_joint_probability_distribution(TS, rang, nbins):
    """
    Assign each node j to a vector of length nbins where each element is the product of its own
    individual_probability_distribution and its neighbors'. P(x) * P(y) <-- as opposed to P(x,y)

    Parameters
    ----------
    TS (np.ndarray): Array consisting of $L$ observations from $N$ sensors.
    rang (list): list of the minimum and maximum value in the time series
    nbins (int): number of bins for the pre-processing step (to yield a discrete probability distribution)

    Returns
    -------
    JointP (dict): a dictionary where the keys are nodes in the graph and the
                         are nbins x nbins arrays corresponding to joint probability vectors
    """

    N, L = TS.shape  # N nodes and L length
    JointP = dict()  # create a dict to put the joint prob distributions into

    for l in range(N):  # for each node
        for j in range(l):  # for each possible edge
            P, _, _ = np.histogram2d(
                TS[j], TS[l], bins=nbins, range=np.array([rang, rang])
            )
            JointP[(j, l)] = P / L

    return JointP 
Example #18
Source File: test_histogram2d.py    From differential-privacy-library with MIT License 5 votes vote down vote up
def test_missing_range(self):
        x = np.array([1, 2, 3, 4, 5])
        y = np.array([5, 7, 1, 5, 9])
        with self.assertWarns(PrivacyLeakWarning):
            res = histogram2d(x, y, epsilon=1, range=[(0, 10), None])
        self.assertIsNotNone(res) 
Example #19
Source File: pid_analysis.py    From flight_review with BSD 3-Clause "New" or "Revised" License 5 votes vote down vote up
def weighted_mode_avr(self, values, weights, vertrange, vertbins):
        ### finds the most common trace and std
        threshold = 0.5  # threshold for std calculation
        filt_width = 7  # width of gaussian smoothing for hist data

        resp_y = np.linspace(vertrange[0], vertrange[-1], vertbins, dtype=np.float64)
        times = np.repeat(np.array([self.time_resp],dtype=np.float64), len(values), axis=0)
        weights = np.repeat(weights, len(values[0]))

        hist2d = np.histogram2d(times.flatten(), values.flatten(),
                                range=[[self.time_resp[0], self.time_resp[-1]], vertrange],
                                bins=[len(times[0]), vertbins], weights=weights.flatten())[0].transpose()
        ### shift outer edges by +-1e-5 (10us) bacause of dtype32. Otherwise different precisions lead to artefacting.
        ### solution to this --> somethings strage here. In outer most edges some bins are doubled, some are empty.
        ### Hence sometimes produces "divide by 0 error" in "/=" operation.

        if hist2d.sum():
            hist2d_sm = gaussian_filter1d(hist2d, filt_width, axis=0, mode='constant')
            hist2d_sm /= np.max(hist2d_sm, 0)


            pixelpos = np.repeat(resp_y.reshape(len(resp_y), 1), len(times[0]), axis=1)
            avr = np.average(pixelpos, 0, weights=hist2d_sm * hist2d_sm)
        else:
            hist2d_sm = hist2d
            avr = np.zeros_like(self.time_resp)
        # only used for monochrome error width
        hist2d[hist2d <= threshold] = 0.
        hist2d[hist2d > threshold] = 0.5 / (vertbins / (vertrange[-1] - vertrange[0]))

        std = np.sum(hist2d, 0)

        return avr, std, [self.time_resp, resp_y, hist2d_sm]

    ### calculates weighted avverage and resulting errors 
Example #20
Source File: test_function_base.py    From lambda-packs with MIT License 5 votes vote down vote up
def test_f32_rounding(self):
        # gh-4799, check that the rounding of the edges works with float32
        x = np.array([276.318359, -69.593948, 21.329449], dtype=np.float32)
        y = np.array([5005.689453, 4481.327637, 6010.369629], dtype=np.float32)
        counts_hist, xedges, yedges = np.histogram2d(x, y, bins=100)
        assert_equal(counts_hist.sum(), 3.) 
Example #21
Source File: plotter.py    From ChainConsumer with MIT License 5 votes vote down vote up
def _get_smoothed_histogram2d(self, chain, param1, param2):  # pragma: no cover
        # No test coverage here because
        smooth = chain.config["smooth"]
        bins = chain.config["bins"]
        kde = chain.config["kde"]

        x = chain.get_data(param1)
        y = chain.get_data(param2)
        w = chain.weights

        if chain.grid:
            binsx = get_grid_bins(x)
            binsy = get_grid_bins(y)
            hist, x_bins, y_bins = np.histogram2d(x, y, bins=[binsx, binsy], weights=w)
        else:
            binsx, smooth = get_smoothed_bins(smooth, bins, x, w, marginalised=False)
            binsy, _ = get_smoothed_bins(smooth, bins, y, w, marginalised=False)
            hist, x_bins, y_bins = np.histogram2d(x, y, bins=[binsx, binsy], weights=w)

        if chain.power is not None:
            hist = hist ** chain.power

        x_centers = 0.5 * (x_bins[:-1] + x_bins[1:])
        y_centers = 0.5 * (y_bins[:-1] + y_bins[1:])

        if kde:
            nn = x_centers.size * 2  # Double samples for KDE because smooth
            x_centers = np.linspace(x_bins.min(), x_bins.max(), nn)
            y_centers = np.linspace(y_bins.min(), y_bins.max(), nn)
            xx, yy = meshgrid(x_centers, y_centers, indexing="ij")
            coords = np.vstack((xx.flatten(), yy.flatten())).T
            data = np.vstack((x, y)).T
            hist = MegKDE(data, w, kde).evaluate(coords).reshape((nn, nn))
            if chain.power is not None:
                hist = hist ** chain.power
        elif smooth:
            hist = gaussian_filter(hist, smooth, mode=self.parent._gauss_mode)

        return hist, x_centers, y_centers 
Example #22
Source File: twod_mjc_env.py    From inverse_rl with MIT License 5 votes vote down vote up
def make_density_map(paths, map_config):
    xs = np.linspace(map_config.xs[0], map_config.xs[1], num=map_config.xres+1)
    ys = np.linspace(map_config.ys[0], map_config.ys[1], num=map_config.yres+1)
    y = paths[:,0]
    x = paths[:,1]
    H, xedges, yedges = np.histogram2d(y, x, bins=(xs, ys))
    H = H.astype(np.float)
    H = H/np.max(H)
    return H.T 
Example #23
Source File: simulate.py    From picasso with MIT License 5 votes vote down vote up
def convertMovie(
    runner,
    photondist,
    structures,
    imagesize,
    frames,
    psf,
    photonrate,
    background,
    noise,
    mode3Dstate,
    cx,
    cy,
):
    edges = range(0, imagesize + 1)

    photonposframe = distphotonsxy(
        runner, photondist, structures, psf, mode3Dstate, cx, cy
    )

    if len(photonposframe) == 0:
        simframe = _np.zeros((imagesize, imagesize))
    else:
        x = photonposframe[:, 0]
        y = photonposframe[:, 1]
        simframe, xedges, yedges = _np.histogram2d(y, x, bins=(edges, edges))
        simframe = _np.flipud(simframe)  # to be consistent with render

    return simframe 
Example #24
Source File: metric.py    From SpaceNet_Off_Nadir_Solutions with Apache License 2.0 5 votes vote down vote up
def calc_score(labels, y_pred):
    if y_pred.sum() == 0 and labels.sum() == 0:
        return 1
    if labels.sum() == 0 and y_pred.sum() > 0 or y_pred.sum() == 0 and labels.sum() > 0:
        return 0

    true_objects = len(np.unique(labels))
    pred_objects = len(np.unique(y_pred))

    intersection = np.histogram2d(labels.flatten(), y_pred.flatten(), bins=(true_objects, pred_objects))[0]

    # Compute areas (needed for finding the union between all objects)
    area_true = np.histogram(labels, bins=true_objects)[0]
    area_pred = np.histogram(y_pred, bins=pred_objects)[0]
    area_true = np.expand_dims(area_true, -1)
    area_pred = np.expand_dims(area_pred, 0)

    # Compute union
    union = area_true + area_pred - intersection

    # Exclude background from the analysis
    intersection = intersection[1:, 1:]
    union = union[1:, 1:]
    union[union == 0] = 1e-9

    # Compute the intersection over union
    iou = intersection / union
    return precision_at(0.5, iou)


# Precision helper function 
Example #25
Source File: metrics.py    From SpaceNet_Off_Nadir_Solutions with Apache License 2.0 5 votes vote down vote up
def calc_iou(gt_masks, predicted_masks, height=768, width=768):
    true_objects = gt_masks.shape[2]
    pred_objects = predicted_masks.shape[2]
    labels = np.zeros((height, width), np.uint16)
    
    for index in range(0, true_objects):
        labels[gt_masks[:, :, index] > 0] = index + 1
    y_true = labels.flatten()
    labels_pred = np.zeros((height, width), np.uint16)
    
    for index in range(0, pred_objects):
        if sum(predicted_masks[:, :, index].shape) == height + width:
            labels_pred[predicted_masks[:, :, index] > 0] = index + 1
    y_pred = labels_pred.flatten()
    
    intersection = np.histogram2d(y_true, y_pred, bins=(true_objects + 1, pred_objects + 1))[0]
    
    area_true = np.histogram(labels, bins=true_objects + 1)[0]
    area_pred = np.histogram(y_pred, bins=pred_objects + 1)[0]
    area_true = np.expand_dims(area_true, -1)
    area_pred = np.expand_dims(area_pred, 0)
    
    # Compute union
    union = area_true + area_pred - intersection
    
    # Exclude background from the analysis
    intersection = intersection[1:, 1:]
    union = union[1:, 1:]
    union[union == 0] = 1e-9
    # print(union)
    # print(intersection)
    iou = intersection / union
    return iou 
Example #26
Source File: test_function_base.py    From ImageFusion with MIT License 5 votes vote down vote up
def test_f32_rounding(self):
        # gh-4799, check that the rounding of the edges works with float32
        x = np.array([276.318359  , -69.593948  , 21.329449], dtype=np.float32)
        y = np.array([5005.689453, 4481.327637, 6010.369629], dtype=np.float32)
        counts_hist, xedges, yedges = np.histogram2d(x, y, bins=100)
        assert_equal(counts_hist.sum(), 3.) 
Example #27
Source File: test_histogram2d.py    From differential-privacy-library with MIT License 5 votes vote down vote up
def test_no_params(self):
        x = np.array([1, 2, 3, 4, 5])
        y = np.array([5, 7, 1, 5, 9])
        with self.assertWarns(PrivacyLeakWarning):
            res = histogram2d(x, y)
        self.assertIsNotNone(res) 
Example #28
Source File: test_stats.py    From scprep with GNU General Public License v3.0 5 votes vote down vote up
def calc_MI(X, Y, bins):
    c_XY = np.histogram2d(X, Y, bins)[0]
    c_X = np.histogram(X, bins)[0]
    c_Y = np.histogram(Y, bins)[0]
    H_X = shan_entropy(c_X)
    H_Y = shan_entropy(c_Y)
    H_XY = shan_entropy(c_XY)
    MI = H_X + H_Y - H_XY
    return MI 
Example #29
Source File: viewsExam.py    From muesli with GNU General Public License v3.0 5 votes vote down vote up
def __call__(self):
        source1 = self.request.GET['source1']
        source2 = self.request.GET['source2']
        data1, max1, name1 = self.getData(source1)
        data2, max2, name2 = self.getData(source2)
        student_ids =  set(data1.keys()).intersection(list(data2.keys()))
        data = np.array([(float(data1[s_id]), float(data2[s_id])) for s_id in student_ids])
        if len(data):
            x = data[:,0]
            y = data[:,1]
        else: x,y = [], []
        bins1 = self.getBins(max1)
        bins2 = self.getBins(max2)

        try:
            hist,xedges,yedges = np.histogram2d(x,y,bins=[bins1, bins2])
        except ValueError:
            # x,y = [] raises ValueErrors in old numpy
            hist = np.array([[0]])
            xedges = [0,1]
            yedges = xedges

        color_stepsize = int(math.ceil(hist.max()/10.0)) or 1
        color_ticks = list(range(0, int(hist.max()) or 1, color_stepsize))
        color_bins = [-color_stepsize/2.0]
        color_bins.extend([t+color_stepsize/2.0 for t in color_ticks])

        norm = matplotlib.colors.BoundaryNorm(color_bins, len(color_ticks))

        extent = [xedges[0], xedges[-1], yedges[0], yedges[-1] ]
        ax = self.fig.add_subplot(111)
        im = ax.imshow(hist.T,extent=extent,interpolation='nearest',
                origin='lower',
                cmap = pyplot.cm.gray_r,
                aspect='auto')
        ax.set_xlabel(name1)
        ax.set_ylabel(name2)
        self.fig.colorbar(im, norm=norm, boundaries=color_bins, ticks=color_ticks)
        corrcoef = np.corrcoef(x, y)[0,1] if len(x)>0 else 0
        ax.set_title("Korrelation = %.2f" % corrcoef)
        return self.createResponse() 
Example #30
Source File: test_histogram2d.py    From differential-privacy-library with MIT License 5 votes vote down vote up
def test_no_range(self):
        x = np.array([1, 2, 3, 4, 5])
        y = np.array([5, 7, 1, 5, 9])
        with self.assertWarns(PrivacyLeakWarning):
            res = histogram2d(x, y, epsilon=1)
        self.assertIsNotNone(res)