diff --git a/docs/notebooks/demo.ipynb b/docs/notebooks/demo.ipynb index 1671a070..306843e6 100644 --- a/docs/notebooks/demo.ipynb +++ b/docs/notebooks/demo.ipynb @@ -90,7 +90,8 @@ }, "outputs": [], "source": [ - "quantiles = P.quantize(percent=10.)\n", + "quants = np.array([0.01,0.02,0.03,0.04,0.05,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,0.95,0.96,0.97,0.98,0.99])\n", + "quantiles = P.quantize(quants=quants)#(percent=10.)\n", "Q = qp.PDF(quantiles=quantiles)" ] }, diff --git a/qp/pdf.py b/qp/pdf.py index 58dfe876..2d1243ab 100644 --- a/qp/pdf.py +++ b/qp/pdf.py @@ -4,7 +4,6 @@ import qp -<<<<<<< HEAD class PDF(object): def __init__(self, truth=None, quantiles=None, histogram=None, @@ -37,8 +36,7 @@ def __init__(self, truth=None, quantiles=None, histogram=None, self.initialized = 'histogram' self.last = self.initialized - if vb and self.truth is None and self.quantiles is None - and self.histogram is None: + if vb and self.truth is None and self.quantiles is None and self.histogram is None: print 'Warning: initializing a PDF object without inputs' self.difs = None self.mids = None @@ -96,12 +94,14 @@ def integrate(self, limits): """ return None - def quantize(self, percent=1., number=None, vb=True): + def quantize(self, quants=None, percent=1., number=None, vb=True): """ Computes an array of evenly-spaced quantiles from the truth. Parameters ---------- + quants: ndarray, float + array of quantile locations as decimals percent: float the separation of the requested quantiles, in percent number: int @@ -123,25 +123,34 @@ def quantize(self, percent=1., number=None, vb=True): object stored in `self.truth`. This calculates the inverse CDF. See `the Scipy docs `_ for details. """ - if number is not None: - # Compute the spacing of the quantiles: - quantum = 1.0 / float(number+1) + if quants is not None: + full_quants = np.append(np.array([0.]),quants) + full_quants = np.append(full_quants,np.array([1.])) + quantdifs = full_quants[1:]-full_quants[:-1] + assert np.sum(quantdifs)==1. + self.quantpoints = quants else: - quantum = percent/100.0 - # Over-write the number of quantiles: - number = np.ceil(100.0 / percent) - 1 - assert number > 0 - - points = np.linspace(0.0+quantum, 1.0-quantum, number) - if vb: print("Calculating quantiles: ", points) + if number is not None: + # Compute the spacing of the quantiles: + quantum = 1.0 / float(number+1) + else: + quantum = percent/100.0 + # Over-write the number of quantiles: + number = np.ceil(100.0 / percent) - 1 + assert number > 0 + + self.quantpoints = np.linspace(0.0+quantum, 1.0-quantum, number) + if vb: print("Calculating quantiles: ", self.quantpoints) if self.truth is not None: - self.quantiles = self.truth.ppf(points) + self.quantiles = self.truth.ppf(self.quantpoints) else: print('New quantiles can only be computed from a truth distribution in this version.') return if vb: print("Result: ", self.quantiles) + self.quantvals = self.evaluate(self.quantiles) self.last = 'quantiles' + self.quantiles = (self.quantiles, self.quantvals) return self.quantiles def histogramize(self, binends=None, nbins=10, binrange=[0., 1.], vb=True): @@ -222,9 +231,10 @@ def interpolate(self, using=None, vb=True): if self.quantiles is None: self.quantiles = self.quantize() - self.difs = self.quantiles[1:]-self.quantiles[:-1] - self.mids = (self.quantiles[1:]+self.quantiles[:-1])/2. - self.vals = (1.0/(len(self.quantiles)+1))/self.difs + #self.difs = self.quantiles[0][1:]-self.quantiles[0][:-1] + #self.mids = (self.quantiles[0][1:]+self.quantiles[0][:-1])/2. + self.mids = self.quantiles[0] + self.vals = self.quantiles[1]#(1.0/(len(self.quantiles)+1))/self.difs if using == 'histogram': # First find the histogram if none exists: @@ -300,7 +310,7 @@ def plot(self, limits, points=None): plt.plot(x, self.truth.pdf(x), color='k', linestyle='-', lw=1.0, alpha=1.0, label='True PDF') if self.quantiles is not None: - plt.vlines(self.quantiles, 0., 1., color='b', linestyle='--', lw=1.0, alpha=1., label='Quantiles') + plt.vlines(self.quantiles[0], np.zeros(len(self.quantiles[1])), self.quantiles[1], color='b', linestyle='--', lw=1.0, alpha=1., label='Quantiles') if points is not None: (grid, qinterpolated) = self.approximate(points, using='quantiles') plt.plot(grid, qinterpolated, color='b', linestyle=':', lw=2.0, alpha=1.0, label='Quantile Interpolated PDF')