I have a few questions about interpolate.InterpolatedUnivariateSpline (IUS).
My ultimate aim, given an IUS object, a lower definite integration limit, and an area, is to obtain the upper definite integration limit by the fastest route possible. The application is to do with probability distribution functions. The IUS represents the PDF. I can quickly work out the CDF using IUS.integral(a, x) = q, where a is the lower limit of support and x is the value for which you want to work out the CDF. However, I would like the quickest way to calculate the percent probability function (PPF). I.e. given q, work out x. One way used in scipy.stats is to use optimize.brentq to find this value - It's known that the CDF is a monotonically increasing function. It strikes me that brentq might not be the fastest way of achieving this. For example, I could cache the integral at each of the datapoints for the IUS, which would immediately allow me to bracket the location of x. Given that I know the degree of the smoothing spline I am hoping to calculate x quickly if I know the polynomial coefficients of the spline in the bracketing interval. This would done by integrating the polynomial and using the roots of the integral. My question is: Is there a way of getting those polynomial coefficients for each interval from the output of the IUS.get_coeffs method? -- _____________________________________
Dr. Andrew Nelson _____________________________________ _______________________________________________ SciPy-User mailing list [hidden email] https://mail.scipy.org/mailman/listinfo/scipy-user |
On Thu, Aug 11, 2016 at 4:40 AM, Andrew Nelson <[hidden email]> wrote:
> I have a few questions about interpolate.InterpolatedUnivariateSpline (IUS). > > My ultimate aim, given an IUS object, a lower definite integration limit, > and an area, is to obtain the upper definite integration limit by the > fastest route possible. > > The application is to do with probability distribution functions. The IUS > represents the PDF. I can quickly work out the CDF using IUS.integral(a, x) > = q, where a is the lower limit of support and x is the value for which you > want to work out the CDF. > > However, I would like the quickest way to calculate the percent probability > function (PPF). I.e. given q, work out x. > > One way used in scipy.stats is to use optimize.brentq to find this value - > It's known that the CDF is a monotonically increasing function. It strikes > me that brentq might not be the fastest way of achieving this. For example, > I could cache the integral at each of the datapoints for the IUS, which > would immediately allow me to bracket the location of x. Given that I know > the degree of the smoothing spline I am hoping to calculate x quickly if I > know the polynomial coefficients of the spline in the bracketing interval. > This would done by integrating the polynomial and using the roots of the > integral. My question is: Is there a way of getting those polynomial > coefficients for each interval from the output of the IUS.get_coeffs method? > First of all, you likely cannot beat an inverse interpolation if raw speed is what you're after (i.e., interpolate cdf(x) vs x). Spefically, to your question about IUS: * note that IUS is *not* guaranteed to be monotonic even if you tabulate what looks like a monotonic function * IUS is defined in the B-spline basis, so .get_coeffs() returns the coeficients of B-spline basis elements, not polynomial representation coefficients. If you want polynomial representation, you can use `PPoly.from_spline()` with knots and coefficients from IUS --- check the docstring of LSQUnivariateSpline for a possible catch with the boundary knots though. If you go that route, you might use `splrep(x, y, s=0)` rather than IUS though, because it's equivalent and gives directly tck (knots, coefficients and order) which `PPoly.from_spline` takes as a parameter. However, if you convert to polynomials, you might as well just use the `CubicSpline` (new in scipy 0.18.0) which constructs the spline directly in the polynomial basis. Unless you want to make sure you get a monotonic interpolator and can tolerate C1 continuity, in which case you can use instead PCHIP or Akima1DInterpolator. For integrals specifically, one thing which might (or might not) help is to use Bernstein polynomials (which e.g. PchipInterpolator is), and use the fact that $\int_0^{1} b_{j, n} = \frac{1}{n+1} $, so that you know the cumulative integrals at each breakpoint. HTH, Evgeni _______________________________________________ SciPy-User mailing list [hidden email] https://mail.scipy.org/mailman/listinfo/scipy-user |
Free forum by Nabble | Edit this page |