pygama.math.functions package#
Statistical distributions for the Pygama pacakge.
Each distribution requires the definition of four vectorized numbafied functions that take an array as an input:
1. nb_dist_pdf(x, shape, mu, sigma)()
Returns the PDF normalized on the support
2. nb_dist_cdf(x, shape, mu, sigma)()
Returns the CDF derived from the PDF that is normalized on the support
3. nb_dist_scaled_pdf(x, area, shape, mu, sigma)()
Returns area*nb_dist_prd
4. nb_dist_scaled_cdf(x, area, shape, mu, sigma)()
Returns area*nb_dist_cdf
NOTE: The order of the arguments of these functions follows the ordering convention from left to right (for whichever are present): x, area, shapes, mu, sigma
Then these four functions are pacakged into a class that subclasses our own pygama_continuous. This distribution class then has 9 required methods.
1. _pdf(x, shape, mu, sigma)()
An overloading of the scipy rv_continuous base class’ method, this is slow, but it enables access to other scipy rv_continuous methods like _rvs
2. _cdf(x, shape, mu, sigma)()
An overloading of the scipy rv_continuous base class’ method, this is slow, but it enables access to other scipy rv_continuous methods like _logcdf
3. get_pdf(x, shape, mu, sigma)()
A direct call to nb_dist_pdf, it is very quick. This PDF is normalized on its support
4. get_cdf(x, shape, mu, sigma)()
A direct call to nb_dist_cdf, it is very quick. This CDF comes from the PDF that is normalized on the support
5. norm_pdf(x, shape, mu, sigma)()
Returns the get_pdf, normalized to unity on the fit range (normalized on [np.amin(x), np.amax(x)]). Needed for unbinned fits
6. norm_cdf(x, shape, mu, sigma)()
Returns get_cdf, normalized on the fit range (normalized on [np.amin(x), np.amax(x)]). Needed for binned fits
7. pdf_ext(x, area, x_lo, x_hi, shape, mu, sigma)()
Returns both the integral of nb_dist_scaled_pdf, as well as nb_dist_scaled_pdf itself. Needed for extended unbinned fits
8. cdf_ext(x, area, shape, mu, sigma)()
Returns a direct call to nb_dist_scaled_cdf itself. Needed for extended binned fits.
9. required_args()
A tuple of the required shape, mu, and sigma parameters
NOTE: the order of the arguments to these functions must follow the convention for whichever are present: x, area, x_lo, x_hi, shape, mu, sigma
pygama_continuous subclasses scipy’s rv_continuous and overloads the class instantiation so that is calls numba_frozen. numba_frozen is a subclass of rv_continuous, and acts exactly like scipy’s rv_continuous_frozen, except that this class also has methods for pygama specific get_pdf, get_cdf, pdf_ext, and cdf_ext. This subclassing is necessary so that we can instantiate frozen distributions that have access to our required pygama class methods.
In addition, the class sum_dists subclasses rv_continuous. This is so that distributions created out of the sum of distributions also have access to scipy rv_continuous methods, such as random sampling. sum_dists also has convenience functions such as draw_pdf, which plots the pdf, and get_req_args which returns a list of the required arguments for the sum of the distribution.
Submodules#
pygama.math.functions.crystal_ball module#
Crystal ball distributions for Pygama
- class pygama.math.functions.crystal_ball.crystal_ball_gen(momtype=1, a=None, b=None, xtol=1e-14, badvalue=None, name=None, longname=None, shapes=None, extradoc=None, seed=None)#
Bases:
pygama_continuous- cdf_ext(x: ndarray, area: float, beta: float, m: float, mu: float, sigma: float) ndarray#
- Return type:
- norm_cdf(x: ndarray, x_lower: float, x_upper: float, beta: float, m: float, mu: float, sigma: float) ndarray#
- Return type:
- norm_pdf(x: ndarray, x_lower: float, x_upper: float, beta: float, m: float, mu: float, sigma: float) ndarray#
- Return type:
- pygama.math.functions.crystal_ball.nb_crystal_ball_cdf(x: ndarray, beta: float, m: float, mu: float, sigma: float) ndarray#
CDF for power-law tail plus gaussian. Its range of support is \(x\in\mathbb{R}, \beta>0, m>1\). It computes:
\[\begin{split}cdf(x, \beta, m, \mu, \sigma)= \begin{cases} NA\sigma\frac{(B-\frac{x-\mu}{\sigma})^{1-m}}{m-1} \quad , \frac{x-\mu}{\sigma} \leq -\beta \\ NA\sigma\frac{(B+\beta)^{1-m}}{m-1} + N\sigma \sqrt{\frac{\pi}{2}}\left(\text{erf}\left(\frac{x-\mu}{\sigma \sqrt{2}}\right)+\text{erf}\left(\frac{\beta}{\sqrt{2}}\right)\right) \quad , \frac{x-\mu}{\sigma} > -\beta \end{cases}\end{split}\]Where
\[\begin{split}A = \frac{m^m}{\beta^m}e^{-\beta^2/2} \\ B = \frac{m}{\beta}-\beta \\ N = \frac{1}{\sigma \frac{m e^{-\beta^2/2}}{\beta(m-1)} + \sigma \sqrt{\frac{\pi}{2}}\left(1+\text{erf}\left(\frac{\beta}{\sqrt{2}}\right)\right)}\end{split}\]As a Numba vectorized function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.crystal_ball.nb_crystal_ball_pdf(x: ndarray, beta: float, m: float, mu: float, sigma: float) ndarray#
PDF of a power-law tail plus gaussian. Its range of support is \(x\in\mathbb{R}, \beta>0, m>1\). It computes:
\[\begin{split}pdf(x, \beta, m, \mu, \sigma) = \begin{cases}NA(B-\frac{x-\mu}{\sigma})^{-m} \quad \frac{x-\mu}{\sigma}\leq -\beta \\ Ne^{-(\frac{x-\mu}{\sigma})^2/2} \quad \frac{x-\mu}{\sigma}>-\beta\end{cases}\end{split}\]Where
\[\begin{split}A = \frac{m^m}{\beta^m}e^{-\beta^2/2} \\ B = \frac{m}{\beta}-\beta \\ N = \frac{1}{\sigma \frac{m e^{-\beta^2/2}}{\beta(m-1)} + \sigma \sqrt{\frac{\pi}{2}}\left(1+\text{erf}\left(\frac{\beta}{\sqrt{2}}\right)\right)}\end{split}\]As a Numba vectorized function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.crystal_ball.nb_crystal_ball_scaled_cdf(x: ndarray, area: float, beta: float, m: float, mu: float, sigma: float) ndarray#
Scaled CDF for power-law tail plus gaussian. Used for extended binned fits. As a Numba vectorized function, it runs slightly faster than ‘out of the box’ functions.
- Parameters:
x (ndarray) – The input data
beta (float) – The point where the cdf changes from power-law to Gaussian
m (float) – The power of the power-law tail
mu (float) – The amount to shift the distribution
sigma (float) – The amount to scale the distribution
area (float) – The number of counts in the distribution
- Return type:
- pygama.math.functions.crystal_ball.nb_crystal_ball_scaled_pdf(x: ndarray, area: float, beta: float, m: float, mu: float, sigma: float) ndarray#
Scaled PDF of a power-law tail plus gaussian. As a Numba vectorized function, it runs slightly faster than ‘out of the box’ functions.
- Parameters:
x (ndarray) – The input data
beta (float) – The point where the pdf changes from power-law to Gaussian
m (float) – The power of the power-law tail
mu (float) – The amount to shift the distribution
sigma (float) – The amount to scale the distribution
area (float) – The number of counts in the distribution
- Return type:
pygama.math.functions.error_function module#
pygama.math.functions.exgauss module#
Exponentially modified Gaussian distributions for pygama
- class pygama.math.functions.exgauss.exgauss_gen(momtype=1, a=None, b=None, xtol=1e-14, badvalue=None, name=None, longname=None, shapes=None, extradoc=None, seed=None)#
Bases:
pygama_continuous- norm_cdf(x: ndarray, x_lower: float, x_upper: float, tau: float, mu: float, sigma: float) ndarray#
- Return type:
- norm_pdf(x: ndarray, x_lower: float, x_upper: float, tau: float, mu: float, sigma: float) ndarray#
- Return type:
- pygama.math.functions.exgauss.nb_exgauss_cdf(x: ndarray, mu: float, sigma: float, tau: float, lower_range: float = inf, upper_range: float = inf) ndarray#
The CDF for a normalized exponentially modified Gaussian. Its range of support is \(x\in(-\infty,\infty)\), \(\tau\in(-\infty,\infty)\) It computes the following CDF:
\[cdf(x,\tau,\mu,\sigma) = \tau\,pdf(x,\tau,\mu,\sigma)+ \frac{\tau}{2|\tau|}\text{erf}\left(\frac{\tau}{|\tau|\sqrt{2}}(\frac{x-\mu}{\sigma})\right) + \frac{1}{2}\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- Parameters:
x (ndarray) – Input data
mu (float) – The centroid of the Gaussian
sigma (float) – The standard deviation of the Gaussian
tau (float) – The characteristic scale of the Gaussian tail
lower_range (float) – Lower bound of the Gaussian with tail
upper_range (float) – Upper bound of the Gaussian with tail
- Return type:
- pygama.math.functions.exgauss.nb_exgauss_pdf(x: ndarray, mu: float, sigma: float, tau: float) ndarray#
Normalized PDF of an exponentially modified Gaussian distribution. Its range of support is \(x\in(-\infty,\infty)\), \(\tau\in(-\infty,\infty)\) Calls either
nb_gauss_tail_exact()ornb_gauss_tail_approx()depending on which is computationally cheaperAs a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- Parameters:
- Return type:
See also
- pygama.math.functions.exgauss.nb_exgauss_scaled_cdf(x: ndarray, mu: float, sigma: float, tau: float, area: float) ndarray#
Scaled CDF of an exponentially modified Gaussian distribution As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.exgauss.nb_exgauss_scaled_pdf(x: ndarray, mu: float, sigma: float, tau: float, area: float) ndarray#
Scaled PDF of an exponentially modified Gaussian distribution Can be used as a component of other fit functions w/args mu,sigma,tau As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.exgauss.nb_gauss_tail_approx(x: ndarray, mu: float, sigma: float, tau: float) ndarray#
Approximate form of a normalized exponentially modified Gaussian PDF As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- Parameters:
- Return type:
See also
- pygama.math.functions.exgauss.nb_gauss_tail_exact(x: float, mu: float, sigma: float, tau: float, tmp: float) float#
Exact form of a normalized exponentially modified Gaussian PDF. It computes the following PDF:
\[pdf(x, \tau,\mu,\sigma) = \frac{1}{2|\tau|}e^{\frac{x-\mu}{\tau}+\frac{\sigma^2}{2\tau^2}}\text{erfc}\left(\frac{\tau(\frac{x-\mu}{\sigma})+\sigma}{|\tau|\sqrt{2}}\right)\]Where \(tmp = \frac{x-\mu}{\tau}+\frac{\sigma^2}{2\tau^2}\) is precomputed in
nb_exgauss_pdf()to save computational time.As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- Parameters:
- Return type:
See also
pygama.math.functions.exponential module#
Exponential distributions for pygama
- class pygama.math.functions.exponential.exponential_gen(momtype=1, a=None, b=None, xtol=1e-14, badvalue=None, name=None, longname=None, shapes=None, extradoc=None, seed=None)#
Bases:
pygama_continuous- norm_cdf(x: ndarray, x_lower: float, x_upper: float, lamb: float, mu: float, sigma: float) ndarray#
- Return type:
- norm_pdf(x: ndarray, x_lower: float, x_upper: float, lamb: float, mu: float, sigma: float) ndarray#
- Return type:
- pygama.math.functions.exponential.nb_exponential_cdf(x: ndarray, lamb: float, mu: float, sigma: float) ndarray#
Normalised exponential cumulative distribution, w/ args: lamb, mu, sigma. Its range of support is \(x\in[0,\infty), \lambda>0\). It computes:
\[\begin{split}cdf(x, \lambda, \mu, \sigma) = \begin{cases} 1-e^{-\lambda\frac{x-\mu}{\sigma}} \quad , \frac{x-\mu}{\sigma} > 0 \\ 0 \quad , \frac{x-\mu}{\sigma}\leq 0 \end{cases}\end{split}\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.exponential.nb_exponential_pdf(x: ndarray, lamb: float, mu: float, sigma: float) ndarray#
Normalised exponential probability density distribution, w/ args: lamb, mu, sigma. Its range of support is \(x\in[0,\infty), \lambda>0\). It computes:
\[\begin{split}pdf(x, \lambda, \mu, \sigma) = \begin{cases} \lambda e^{-\lambda\frac{x-\mu}{\sigma}} \quad , \frac{x-\mu}{\sigma}\geq 0 \\ 0 \quad , \frac{x-\mu}{\sigma}<0 \end{cases}\end{split}\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.exponential.nb_exponential_scaled_cdf(x: ndarray, lamb: float, mu: float, sigma: float, area: float) ndarray#
Exponential cdf scaled by the area, used for extended binned fits As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
pygama.math.functions.gauss module#
Gaussian distributions for pygama
- class pygama.math.functions.gauss.gaussian_gen(momtype=1, a=None, b=None, xtol=1e-14, badvalue=None, name=None, longname=None, shapes=None, extradoc=None, seed=None)#
Bases:
pygama_continuous
- pygama.math.functions.gauss.nb_gauss(x: ndarray, mu: float, sigma: float) ndarray#
Gaussian, unnormalised for use in building PDFs As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.gauss.nb_gauss_amp(x: ndarray, mu: float, sigma: float, a: float) ndarray#
Gaussian with height as a parameter for FWHM etc. As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.gauss.nb_gauss_cdf(x: ndarray, mu: float, sigma: float) ndarray#
Gaussian CDF, w/ args: mu, sigma. The support is \((-\infty, \infty)\)
\[cdf(x, \mu,\sigma) = \frac{1}{2}\left[1+\text{erf}(\frac{x-\mu}{\sigma\sqrt{2}})\right]\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.gauss.nb_gauss_pdf(x: ndarray, mu: float, sigma: float) ndarray#
Normalised Gaussian PDF, w/ args: mu, sigma. The support is \((-\infty, \infty)\)
\[pdf(x, \mu, \sigma) = \frac{1}{\sqrt{2\pi}}e^{(\frac{x-\mu}{\sigma}^2)/2}\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.gauss.nb_gauss_scaled_cdf(x: ndarray, mu: float, sigma: float, area: float) ndarray#
Gaussian CDF scaled by the number of signal counts for extended binned fits As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
pygama.math.functions.gauss_on_exgauss module#
- class pygama.math.functions.gauss_on_exgauss.gauss_on_exgauss_gen#
Bases:
sum_distsProvide a convenience function for a Gaussian on top of an extended Gaussian. The correct usage is
gauss_on_exgauss.pdf(x, pars = [mu, sigma, htail, tau])()- Parameters:
mu – The location and scale of the first Gaussian
htail – The height of the tail
tau – The characteristic scale of the extended Gaussian tail
- Returns:
gauss_on_exgauss – A subclass of
sum_distsandrv_continuous, has methods ofpdf(),cdf(), etc.
Notes
The extended Gaussian distribution shares the mu, sigma with the Gaussian
- _link_pars(shape_par_idx, area_idx, frac_idx, total_area_idx, params, areas, fracs, total_area)#
pygama.math.functions.gauss_on_linear module#
- class pygama.math.functions.gauss_on_linear.gauss_on_linear_gen#
Bases:
sum_distsProvide a convenience function for a Gaussian on top of a linear distribution.
- Parameters:
area1 – The area of the Gaussian distribution
mu – The location and scale of the first Gaussian
sigma – The location and scale of the first Gaussian
area2 – The area of the linear distributions respectively
m – The slope of the linear distribution
b – The y intercept of the linear distribution
- Returns:
gauss_on_linear – A subclass of sum_dists and rv_continuous, has methods of pdf, cdf, etc.
Notes
link_pars automatically makes sure that the linear distribution is only evaluated on the range of the x array eventually passed Parameters must be in the order (area1, mu, sigma, area2, m, b)
- _link_pars(shape_par_idx, area_idx, frac_idx, total_area_idx, params, areas, fracs, total_area)#
pygama.math.functions.gauss_on_step module#
- class pygama.math.functions.gauss_on_step.gauss_on_step_gen#
Bases:
sum_distsProvide a convenience function for a Gaussian on top of a step distribution.
- Parameters:
area1 – The area of the Gaussian distribution
mu – The location and scale of the first Gaussian, and step function
sigma – The location and scale of the first Gaussian, and step function
area2 – The area of the step function
hstep – The height of the step function
lower_range – The lower and upper bounds of the support of the step function
upper_range – The lower and upper bounds of the support of the step function
- Returns:
gauss_on_step – A subclass of sum_dists and rv_continuous, has methods of pdf, cdf, etc.
Notes
The step function shares a mu and sigma with the step function The parameter array must have ordering (area1, mu, sigma, area2, hstep, lower_range, upper_range)
- _cdf(x, params)#
- get_cdf(x, params)#
pygama.math.functions.gauss_on_uniform module#
- class pygama.math.functions.gauss_on_uniform.gauss_on_uniform_gen#
Bases:
sum_distsProvide a convenience function for a Gaussian on top of a uniform distribution.
- Parameters:
area1 – The area of the Gaussian distribution
mu – The location and scale of the first Gaussian
sigma – The location and scale of the first Gaussian
area2 – The area of the uniform distributions respectively
- Returns:
gauss_on_uniform – A subclass of sum_dists and rv_continuous, has methods of pdf, cdf, etc.
Notes
link_pars automatically makes sure that the uniform distribution is only evaluated on the range of the x array eventually passed
- _link_pars(shape_par_idx, area_idx, frac_idx, total_area_idx, params, areas, fracs, total_area)#
pygama.math.functions.hpge_peak module#
- class pygama.math.functions.hpge_peak.hpge_peak_gen#
Bases:
sum_distsProvide a convenience function for the HPGe peak shape.
A HPGe peak consists of a Gaussian on an Exgauss on a step function.
\[PDF = n_sig*((1-htail)*gauss + htail*exgauss) + n_bkg*step\]Called with
hpge_peak.get_pdf(x, params=[n_sig, mu, sigma, htail, tau, n_bkg, hstep, lower_range, upper_range])
- Parameters:
n_sig – The area of the gauss on exgauss
mu – The centroid of the Gaussian
sigma – The standard deviation of the Gaussian
htail – The height of the Gaussian tail
tau – The characteristic scale of the Gaussian tail
n_bkg – The area of the step background
hstep – The height of the step function background
lower_range – Lower bound of the step function
upper_range – Upper bound of the step function
- Returns:
hpge_peak – A subclass of sum_dists and rv_continuous, has methods of pdf, cdf, etc.
Notes
The extended Gaussian distribution and the step distribution share the mu, sigma with the Gaussian
- _link_pars(shape_par_idx, area_idx, frac_idx, total_area_idx, params, areas, fracs, total_area)#
pygama.math.functions.linear module#
- class pygama.math.functions.linear.linear_gen(momtype=1, a=None, b=None, xtol=1e-14, badvalue=None, name=None, longname=None, shapes=None, extradoc=None, seed=None)#
Bases:
pygama_continuous- cdf_ext(x: ndarray, area: float, x_lower: float, x_upper: float, m: float, b: float) ndarray#
- Return type:
- pygama.math.functions.linear.nb_linear_cdf(x: ndarray, x_lower: float, x_upper: float, m: float, b: float) ndarray#
Normalised linear cumulative density function, w/ args: m, b. Its range of support is \(x\in(x_{lower},x_{upper})\). If \(x_{lower} = np.inf\) and \(x_{upper} = np.inf\), then the function takes \(x_{upper} = :func:`np.min(x)\) and \(x_{upper} = :func:`np.amax(x)\) It computes:
\[cdf(x, x_{lower}, x_{upper}, m, b) = \frac{\frac{m}{2}(x^2-x_{lower}^2)+b(x-x_{lower})}{\frac{m}{2}(x_{upper}^2-x_{lower}^2)+b(x_{upper}-x_{lower})}\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.linear.nb_linear_pdf(x: ndarray, x_lower: float, x_upper: float, m: float, b: float) ndarray#
Normalised linear probability density function, w/ args: m, b. Its range of support is \(x\in(x_{lower},x_{upper})\). If \(x_{lower} = np.inf\) and \(x_{upper} = np.inf\), then the function takes \(x_{upper} = :func:`np.min(x)\) and \(x_{upper} = :func:`np.amax(x)\) It computes:
\[pdf(x, x_{lower}, x_{upper}, m, b) = \frac{mx+b}{\frac{m}{2}(x_{upper}^2-x_{lower}^2)+b(x_{upper}-x_{lower})}\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.linear.nb_linear_scaled_cdf(x: ndarray, x_lower: float, x_upper: float, m: float, b: float, area: float) ndarray#
Linear cdf scaled by the area for extended binned fits As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- Parameters:
- Return type:
pygama.math.functions.moyal module#
Moyal distributions for pygama
- class pygama.math.functions.moyal.moyal_gen(momtype=1, a=None, b=None, xtol=1e-14, badvalue=None, name=None, longname=None, shapes=None, extradoc=None, seed=None)#
Bases:
pygama_continuous
- pygama.math.functions.moyal.nb_moyal_cdf(x: ndarray, mu: float, sigma: float) ndarray#
Normalised Moyal cumulative distribution, w/ args: mu, sigma. Its support is \(x\in\mathbb{R}\) It computes:
\[cdf(x, \mu, \sigma) = \text{erfc}\left(\exp\left(-\frac{x-\mu}{2\sigma}\right)/\sqrt{2}\right)\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.moyal.nb_moyal_pdf(x: ndarray, mu: float, sigma: float) ndarray#
Normalised Moyal probability distribution function, w/ args: mu, sigma. Its support is \(x\in\mathbb{R}\) It computes:
\[pdf(x, \mu, \sigma) = \frac{\exp\left(-\left(\frac{x-\mu}{\sigma}+\exp\left(-\frac{x-\mu}{\sigma}\right)\right)/2\right)}{\sigma\sqrt{2\pi}}\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.moyal.nb_moyal_scaled_cdf(x: ndarray, mu: float, sigma: float, area: float) ndarray#
Moyal cdf scaled by the area, used for extended binned fits As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
pygama.math.functions.poisson module#
Poisson distributions for pygama
- pygama.math.functions.poisson.factorial(nn)#
- pygama.math.functions.poisson.nb_poisson_cdf(x: ndarray, lamb: float, mu: int) ndarray#
Normalised Poisson cumulative distribution, w/ args: lamb, mu. As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
\[cdf(x, \lambda, \mu) = e^{-\lambda}\sum_{j=0}^{\lfloor x-\mu \rfloor}\frac{\lambda^j}{j!}\]
- pygama.math.functions.poisson.nb_poisson_pmf(x: ndarray, lamb: float, mu: int) ndarray#
Normalised Poisson distribution, w/ args: lamb, mu. The range of support is \(\mathbb{N}\), with \(lamb\) \(\in (0,\infty)\), \(\mu \in \mathbb{N}\) As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
\[pmf(x, \lambda, \mu) = \frac{\lambda^{x-\mu} e^{-\lambda}}{(x-\mu)!}\]
- pygama.math.functions.poisson.nb_poisson_scaled_cdf(x: ndarray, lamb: float, mu: int, area: float) ndarray#
Poisson cdf scaled by the number of signal counts for extended binned fits As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.poisson.nb_poisson_scaled_pmf(x: ndarray, lamb: float, mu: int, area: float) ndarray#
Scaled Poisson probability distribution, w/ args: lamb, mu. As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- class pygama.math.functions.poisson.poisson_gen(a=0, b=inf, name=None, badvalue=None, moment_tol=1e-08, values=None, inc=1, longname=None, shapes=None, extradoc=None, seed=None)#
Bases:
rv_discrete- _cdf(x, lamb)#
- _pmf(x, lamb)#
- cdf_ext(x, area, lamb, mu)#
- get_cdf(x, lamb, mu)#
- get_pmf(x, lamb, mu)#
- pmf_ext(x, area, x_lo, x_hi, lamb, mu)#
pygama.math.functions.polynomial module#
- pygama.math.functions.polynomial.nb_poly(x: ndarray, pars: ndarray) ndarray#
A polynomial function with pars following the polyfit convention. It computes:
\[y = a_n x^n + ... + a_1 x + a_0\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- Parameters:
- Returns:
result – The polynomial defined by pars at the evaluated at x
- Return type:
Notes
This follows the
numpy.polyfit()convention of \(a_n x^n + ... + a_1 x + a_0\)
pygama.math.functions.pygama_continuous module#
Subclasses of scipy.stat’s rv_continuous, allows for freezing of pygama numba-fied distributions for faster functions Example:
>>> moyal = moyal_gen(name='moyal')
>>> moyal = moyal(1, 2)
>>> moyal.get_pdf([-1,2,3]) # a direct call to the faster numbafied method
>>> moyal.pdf([-1,2,3]) # a call to the slower scipy method, allows access to other scipy methods like .rvs
>>> moyal.rvs(100) # Can access scipy methods!
NOTE: the dist.pdf method is a slow scipy method, and the dist.get_pdf method is fast
- class pygama.math.functions.pygama_continuous.numba_frozen(dist, *args, **kwds)#
Bases:
rv_frozenEssentially, an overloading of the definition of rv_continuous_frozen so that our pygama class instantiations have access to both the slower scipy methods, as well as the faster pygama methods (like get_pdf)
- cdf_ext(x: ndarray, area: float) ndarray#
Direct access to numba-fied extended pdfs, a fast function
- Returns:
scaled support normalized cdf
- Return type:
- get_cdf(x: ndarray) ndarray#
Direct access to numba-fied cdfs, a fast function
- Returns:
support normalized cdf
- Return type:
- get_pdf(x: ndarray) ndarray#
Direct access to numba-fied pdfs, a fast function
- Returns:
support normalized pdf
- Return type:
- logpdf(x: ndarray) ndarray#
Overloading of scipy’s logpdf function, this is a slow function
- Return type:
- norm_cdf(x, x_lower, x_upper)#
Direct access to numba-fied cdfs, a fast function. Normalized on a fit range
- Returns:
fit-range normalized cdf
- norm_pdf(x, x_lower, x_upper)#
Direct access to numba-fied pdfs, a fast function. Normalized on a fit range
- Returns:
fit-range normalized pdf
- class pygama.math.functions.pygama_continuous.pygama_continuous(momtype=1, a=None, b=None, xtol=1e-14, badvalue=None, name=None, longname=None, shapes=None, extradoc=None, seed=None)#
Bases:
rv_continuousSubclass rv_continuous, and modify the instantiation so that we call an overloaded version of rv_continuous_frozen that has direct access to pygama numbafied functions
- _norm_cdf(x, x_lower, x_upper, *args, **kwds)#
Derive a cdf from a pdf that is normalized on a subset of its support, typically over a fit-range.
- Parameters:
x – The data
x_lower – The lower range to normalize on
x_upper – The upper range to normalize on
args – The shape and location and scale parameters of a specific distribution
- Returns:
norm_pdf – The pdf that is normalized on a smaller range
Notes
If upper_range and lower_range are both
np.inf(), then the function automatically takes x_lower =:func:np.amin(x), x_upper=:func:np.amax(x) We also need to overload this in every subclass because we want functions to be able to introspect the shape and location/scale names. For distributions that are only defined on a limited range, like with lower_range, upper_range, we don’t need to call these, instead just call the normal pdf.
- _norm_pdf(x, x_lower, x_upper, *args, **kwds)#
Normalize a pdf on a subset of its support, typically over a fit-range.
- Parameters:
x – The data
x_lower – The lower range to normalize on
x_upper – The upper range to normalize on
args – The shape and location and scale parameters of a specific distribution
- Returns:
norm_pdf – The pdf that is normalized on a smaller range
Notes
If upper_range and lower_range are both
np.inf(), then the function automatically takes x_lower =:func:np.amin(x), x_upper=:func:np.amax(x) We also need to overload this in every subclass because we want functions to be able to introspect the shape and location/scale names For distributions that are only defined on a limited range, like with lower_range, upper_range, we don’t need to call these, instead just call the normal pdf.
pygama.math.functions.step module#
Step distributions for pygama
- pygama.math.functions.step.nb_step_cdf(x: ndarray, mu: float, sigma: float, hstep: float, lower_range: float = inf, upper_range: float = inf) ndarray#
Normalized CDF for step function w/args mu, sigma, hstep, lower_range, upper_range. Its range of support is \(x\in\) (lower_range, upper_range). It computes:
\[cdf(x, hstep, \mu, \sigma, \text{lower_range}, \text{upper_range}) = cdf(y=\frac{x-\mu}{\sigma}, hstep, \text{lower_range}, \text{upper_range}) = \frac{(y-y_{min}) +hstep\left(y\text{erf}(\frac{y}{\sqrt{2}})+\sqrt{\frac{2}{\pi}}e^{-y^2/2}-y_{min}\text{erf}(\frac{y_{min}}{\sqrt{2}})+\sqrt{\frac{2}{\pi}}e^{-y_{min}^2/2}\right)}{\sigma\left[(y_{max}-y_{min}) +hstep\left(y_{max}\text{erf}(\frac{y_{max}}{\sqrt{2}})+\sqrt{\frac{2}{\pi}}e^{-y_{max}^2/2}-y_{min}\text{erf}(\frac{y_{min}}{\sqrt{2}})+\sqrt{\frac{2}{\pi}}e^{-y_{min}^2/2}\right)\right] }\]Where \(y_{max} = \frac{\text{upper_range} - \mu}{\sigma}, y_{min} = \frac{\text{lower_range} - \mu}{\sigma}\). As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- Parameters:
x (ndarray) – The input data
mu (float) – The location of the step
sigma (float) – The “width” of the step, because we are using an error function to define it
hstep (float) – The height of the step
lower_range (float) – The lower range on which to normalize the step PDF, default is to normalize from min to max x values
upper_range (float) – The upper range on which to normalize the step PDF
- Return type:
- pygama.math.functions.step.nb_step_int(x: ndarray, mu: float, sigma: float, hstep: float) ndarray#
Integral of step function w/args mu, sigma, hstep. It computes:
\[\int pdf(x, hstep, \mu, \sigma)\, dx = \sigma\left(\frac{x-\mu}{\sigma} + hstep \left(\frac{x-\mu}{\sigma}\text{erf}\left(\frac{x-\mu}{\sigma\sqrt{2}}\right) + \sqrt{\frac{2}{\pi}}\exp\left(-(\frac{x-\mu}{\sigma})^2/2\right) \right)\right)\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.step.nb_step_pdf(x: ndarray, mu: float, sigma: float, hstep: float, lower_range: float = inf, upper_range: float = inf) ndarray#
Normalised step function w/args mu, sigma, hstep, lower_range, upper_range. Its range of support is \(x\in\) (lower_range, upper_range). It computes:
\[pdf(x, hstep, \mu, \sigma , \text{lower_range}, \text{upper_range}) = pdf(y=\frac{x-\mu}{\sigma}, step, \text{lower_range}, \text{upper_range}) = \frac{1+hstep\text{erf}\left(\frac{x-\mu}{\sigma\sqrt{2}}\right)}{\sigma\left[(y-y_{min}) +hstep\left(y\text{erf}(\frac{y}{\sqrt{2}})+\sqrt{\frac{2}{\pi}}e^{-y^2/2}-y_{min}\text{erf}(\frac{y_{min}}{\sqrt{2}})+\sqrt{\frac{2}{\pi}}e^{-y_{min}^2/2}\right)\right]}\]Where \(y_{max} = \frac{\text{upper_range} - \mu}{\sigma}, y_{min} = \frac{\text{lower_range} - \mu}{\sigma}\). As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- Parameters:
x (ndarray) – The input data
mu (float) – The location of the step
sigma (float) – The “width” of the step, because we are using an error function to define it
hstep (float) – The height of the step
lower_range (float) – The lower range on which to normalize the step PDF, default is to normalize from min to max x values
upper_range (float) – The upper range on which to normalize the step PDF
- Return type:
- pygama.math.functions.step.nb_step_scaled_cdf(x: ndarray, mu: float, sigma: float, hstep: float, lower_range: float, upper_range: float, area: float) ndarray#
Scaled step function CDF w/args mu, sigma, hstep, lower_range, upper_range, area. Used for extended binned fits. As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- Parameters:
x (ndarray) – The input data
mu (float) – The location of the step
sigma (float) – The “width” of the step, because we are using an error function to define it
hstep (float) – The height of the step
lower_range (float) – The lower range on which to normalize the step PDF, default is to normalize from min to max x values
upper_range (float) – The upper range on which to normalize the step PDF
area (float) – The number to scale the distribution by
- Return type:
- pygama.math.functions.step.nb_step_scaled_pdf(x: ndarray, mu: float, sigma: float, hstep: float, lower_range: float, upper_range: float, area: float) ndarray#
Scaled step function pdf w/args mu, sigma, hstep, lower_range, upper_range, area. As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- Parameters:
x (ndarray) – The input data
mu (float) – The location of the step
sigma (float) – The “width” of the step, because we are using an error function to define it
hstep (float) – The height of the step
lower_range (float) – The lower range on which to normalize the step PDF, default is to normalize from min to max x values
upper_range (float) – The upper range on which to normalize the step PDF
area (float) – The number to scale the distribution by
- Return type:
- pygama.math.functions.step.nb_unnorm_step_pdf(x: float, mu: float, sigma: float, hstep: float) float#
Unnormalised step function for use in pdfs. It computes:
\[pdf(x, hstep, \mu, \sigma) = 1+ hstep\text{erf}\left(\frac{x-\mu}{\sigma\sqrt{2}}\right)\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- class pygama.math.functions.step.step_gen(momtype=1, a=None, b=None, xtol=1e-14, badvalue=None, name=None, longname=None, shapes=None, extradoc=None, seed=None)#
Bases:
pygama_continuous- cdf_ext(x: ndarray, area: float, hstep: float, lower_range: float = inf, upper_range: float = inf, mu: float = 0, sigma: float = 1) ndarray#
- Return type:
- get_cdf(x: ndarray, hstep: float, lower_range: float = inf, upper_range: float = inf, mu: float = 0, sigma: float = 1) ndarray#
- Return type:
- get_pdf(x: ndarray, hstep: float, lower_range: float = inf, upper_range: float = inf, mu: float = 0, sigma: float = 1) ndarray#
- Return type:
- norm_cdf(x: ndarray, x_lower: float, x_upper: float, hstep: float, mu: float, sigma: float) ndarray#
- Return type:
- norm_pdf(x: ndarray, x_lower: float, x_upper: float, hstep: float, mu: float, sigma: float) ndarray#
- Return type:
pygama.math.functions.sum_dists module#
A class that creates the sum of distributions, with methods for scipy computed pdfs() and cdfs(), as well as fast get_pdfs(), and pdf_ext()
mu1, mu2, sigma = range(3)
moyal_add = sum_dists(moyal, [mu1, sigma], moyal, [mu2, sigma], areas = [0.75, 0.25]) # create two moyals that share a sigma,
x = np.arange(-10,10)
pars = np.array([1,2,2]) # corresponds to mu1 = 1, mu2 = 2, sigma = 2
moyal_add.pdf(x, pars)
moyal_add.draw_pdf(x, pars)
moyal_add.get_req_args()
- pygama.math.functions.sum_dists._precompute_fracs(fracs, dists)#
When initializing a summed dist, check the number of fracs provided. If len(fracs) = len(dists)-1, then the remaining frac must be 1-sum(fracs) If len(fracs) = len(dists), return the fracs as given Otherwise, raise an error.
- pygama.math.functions.sum_dists._precompute_shape_par_idx(shape_par_idx, dists)#
Check that each distribution has received the correct number of shape and location parameters
- Parameters:
shape_par_idx – An array of arrays, each containing the indices of the eventual parameters the function call will take from a single parameter array
dists – An array of pygama_continuous distributions, each must have the required_args method!
- Returns:
shape_par_idx – Returned only if every distribution has the correct number of required arguments
- pygama.math.functions.sum_dists.get_idx(idx_array, frac_flag)#
Create separate arrays of parameter indices for the shape parameters, the area, the fracs, and the total area
This is done because the user can pass a combination of areas and fracs during a method call, or can send them in during class initialization. So we need to know what array the actual values of the fracs and areas will be in.
This will return None for an array index whose values are not located in the parameter array passed to a method call This allows
link_pars()to use either the default 1s for the missing areas/fracs/total_area, or user-specified values (passed using keywords).
- class pygama.math.functions.sum_dists.sum_dist_frozen(dist, *args, **kwds)#
Bases:
rv_frozenEssentially, an overloading of the definition of rv_continuous_frozen so that our pygama class instantiations have access to both the slower scipy methods, as well as the faster pygama methods (like get_pdf)
- cdf_ext(x)#
- get_cdf(x)#
- get_pdf(x)#
- get_req_args()#
- logpdf(x)#
- norm_cdf(x)#
- norm_pdf(x)#
- pdf(x)#
- pdf_ext(x)#
- class pygama.math.functions.sum_dists.sum_dists(*args, **kwargs)#
Bases:
rv_continuousInitialize and rv_continuous method so that we gain access to computable methods. Also precompute the fracs of the linear combination, as well as the support of the sum of the distributions.
The correct way to initialize is sum_dists(d1, p1, d2, p2, …, dn, pn, frac_flag) Where d_i is a distribution and p_i is a parameter index array for that distribution.
Parameter index arrays contain indices that slice a single parameter array that is passed to method calls. For example, if the user will eventually pass parameters=[mu, sigma, frac1, frac2] to function.get_pdf(x, parameters) and the first distribution takes mu, sigma, frac1 as its parameters, then p1=[0,1,2]. If the second distribution takes mu, sigma, frac2 then its parameter index array would be p2=[0, 1, 3] because frac2 is the index 3 entry in parameters.
Each par array can contain [shape, mu, sigma, area, frac] with area and frac being optional, and depend on the flag sent to the constructor, but must be placed in that order.
There are 4 flag options:
flag = “areas”, areas passed as fit variables, global fracs optional (but can be passed by kwarg, default all to 1)
- flag = “fracs”, fracs and total area passed as fit variables, global areas optional (but can be passed by kwarg, default all to 1)
by default, the total area is just a parameter of the last distribution. I.e. dist1, [shapes, mu, sigma, frac1], …, distn, [shapes, mu, sigma, frac2, total_area]
flag = “both”, fracs and areas passed as fit variables. dist1, [shapes, mu, sigma, area1, frac1], …, distn, [shapes, mu, sigma, area_n, frac_n, total_area]
flag = None, no areas or fracs passed as fit variables, all default to 1
Notes
dists are the unfrozen pygama distributions
- _argcheck(*args)#
Default check for correct values on args and keywords. Returns condition array of 1’s where arguments are correct and 0’s where they are not.
Pygama needs to overload this to extend it to arrays of arguments, based on how the sum_dists class is built
- _argcheck_rvs(*args, **kwargs)#
- _attach_methods()#
Attaches dynamically created methods to the rv_continuous instance.
- _cdf(x, params)#
- _construct_argparser(meths_to_inspect, locscale_in, locscale_out)#
Construct the parser string for the shape arguments. This method should be called in __init__ of a class for each distribution. It creates the _parse_arg_template attribute that is then used by _attach_argparser_methods to dynamically create and attach the _parse_args, _parse_args_stats, _parse_args_rvs methods to the instance. If self.shapes is a non-empty string, interprets it as a comma-separated list of shape parameters. Otherwise inspects the call signatures of meths_to_inspect and constructs the argument-parsing functions from these. In this case also sets shapes and numargs.
- _link_pars(shape_par_idx, area_idx, frac_idx, total_area_idx, params, areas, fracs, total_area)#
- _pdf(x, params)#
- _ppf(q, shift, scale, *args)#
- _ppf_to_solve(x, q, shift, scale, *args)#
- _rvs(*args, size=None, random_state=None, low=0.0, high=1.0)#
- cdf_ext(x, params)#
- get_cdf(x, params)#
- get_fwhm(pars: ndarray, cov: Optional[ndarray] = None) tuple#
Get the fwhm value from the output of a fit quickly Need to overload this to use hpge_peak_fwhm (to avoid a circular import) for when self is an hpge peak, and otherwise returns 2sqrt(2log(2))*sigma
- get_mu(pars: ndarray, cov: Optional[ndarray] = None, errors: Optional[ndarray] = None) tuple#
Get the mu value from the output of a fit quickly
- get_pdf(x, params)#
- get_req_args()#
This is a default function to get the required args from the distributions in the instance of
sum_dists. This should be overloaded in user defined functions to ensure thatget_mu(),get_fwhm(),get_total_events()
- get_total_events(pars: ndarray, cov: Optional[ndarray] = None, errors: Optional[ndarray] = None) tuple#
Get the total events from the output of an extended fit quickly The total number of events come from the sum of the area of the signal and the area of the background components
- link_pars(shape_par_idx, area_idx, frac_idx, total_area_idx, params, areas, fracs, total_area)#
Overload this if you want to link specific parameters together! Need to overload by first calling pars, areas, fracs, total_area = super()._link_pars
- def _link_pars(self, shape_par_idx, area_idx, frac_idx, total_area_idx, params, areas, fracs, total_area):
shape_pars, cum_len, areas, fracs, total_area = super()._link_pars(shape_par_idx, area_idx, frac_idx, total_area_idx, params, areas, fracs, total_area) fracs[1]= 1 - frac1 total_area[0] = 1 return shape_pars, cum_len, areas, fracs, total_area (mu, sigma, area1, frac1, tau, area2) = range(6)
gauss_on_exgauss = sum_dists(d1 = gauss, p1 = [mu, sigma, area1, frac1], d2 = exgauss, p2 = [tau, mu, sigma, frac1, area2]) pars = [1,2,10,0.75, 0.1, 30] x = np.arange(-10,10) gauss_on_exgauss.pdf(x, pars)
This computes \(10*0.75*gauss(x, mu=1, sigma=2) + 30*(1-0.75)*exgauss(tau = 0.1, mu=1, sigma=2)\) Which is what we want it to do…
- norm_cdf(x, *args, **kwds)#
Derive a cdf from a pdf that is normalized on a subset of its support, typically over a fit-range.
- Parameters:
x – The data
args – The parameter array, can contain x_lo and x_hi at the start of the array
- Returns:
norm_pdf – The pdf that is normalized on a smaller range
Notes
If upper_range and lower_range are both
np.inf(), then the function automatically takes x_lower =:func:np.amin(x), x_upper=:func:np.amax(x) We also need to overload this in every subclass because we want functions to be able to introspect the shape and location/scale names. For distributions that are only defined on a limited range, like with lower_range, upper_range, we don’t need to call these, instead just call the normal pdf.
- norm_pdf(x, *args, **kwds)#
Normalize a pdf on a subset of its support, typically over a fit-range.
- Parameters:
x – The data
args – The parameter array, can start with x_lo, x_hi
- Returns:
norm_pdf – The pdf that is normalized on a smaller range
Notes
If upper_range and lower_range are both
np.inf(), then the function automatically takes x_lower =:func:np.amin(x), x_upper=:func:np.amax(x) We also need to overload this in every subclass because we want functions to be able to introspect the shape and location/scale names For distributions that are only defined on a limited range, like with lower_range, upper_range, we don’t need to call these, instead just call the normal pdf.
- pdf_ext(x, params)#
Extended probability distribution.
Notes
Params can have x_lo, x_hi as the first two parameters, and the rest of the shape parameters in the rest of the array If params is the same length as the required params, the x_lo and x_hi are interpreted to be the max and min of the input array x
- rvs(*args, **kwds)#
Random variates of given type.
- Parameters:
arg1 (array_like) – The shape parameter(s) for the distribution (see docstring of the instance object for more information).
arg2 (array_like) – The shape parameter(s) for the distribution (see docstring of the instance object for more information).
arg3 (array_like) – The shape parameter(s) for the distribution (see docstring of the instance object for more information).
... (array_like) – The shape parameter(s) for the distribution (see docstring of the instance object for more information).
loc (array_like, optional) – Location parameter (default=0).
scale (array_like, optional) – Scale parameter (default=1).
size (int or tuple of ints, optional) – Defining number of random variates (default is 1).
low (float) – The lower bound to sample a uniform distribution
high (float) – The upper bound to sample a uniform distribution
random_state ({None, int, numpy.random.Generator,) – numpy.random.RandomState}, optional If seed is None (or np.random), the numpy.random.RandomState singleton is used. If seed is an int, a new
RandomStateinstance is used, seeded with seed. If seed is already aGeneratororRandomStateinstance then that instance is used.
- Returns:
rvs (ndarray or scalar) – Random variates of given size.
pygama.math.functions.triple_gauss_on_double_step module#
- class pygama.math.functions.triple_gauss_on_double_step.triple_gauss_on_double_step_gen#
Bases:
sum_distsProvide a convenience function for three gaussians on two steps
- Parameters:
area1 – The area, location, and scale of the first Gaussian
mu_1 – The area, location, and scale of the first Gaussian
sigma_1 – The area, location, and scale of the first Gaussian
area2 – The area, location, and scale of the second Gaussian
mu_2 – The area, location, and scale of the second Gaussian
sigma_2 – The area, location, and scale of the second Gaussian
area3 – The area, location, and scale of the third Gaussian
mu_3 – The area, location, and scale of the third Gaussian
sigma_3 – The area, location, and scale of the third Gaussian
area4 – The area and height of the first step function
hstep_1 – The area and height of the first step function
area5 – The area and height of the second step function
hstep_2 – The area and height of the second step function
lower_range – The lower range to compute the normalization of the step functions
upper_range – The upper range to compute the normalization of the step functions
Example
triple_gauss_on_double_step.get_pdf(x, pars = [n_sig1, mu1, sigma1, n_sig2, mu2, sigma2, n_sig3, mu3, sigma3, n_bkg1, hstep1, n_bkg2, hstep2, lower_range, upper_range])
- Returns:
triple_gauss_on_double_step – A subclass of sum_dists and rv_continuous, has methods of pdf, cdf, etc.
Notes
The first step function shares the mu_1, sigma_1 with the first Gaussian, and the second step function shares the mu_2, sigma_2 with the second Gaussian
pygama.math.functions.uniform module#
Uniform distributions for pygama
- pygama.math.functions.uniform.nb_uniform_cdf(x: ndarray, a: float, b: float) ndarray#
Normalised uniform cumulative distribution, w/ args: a, b. Its range of support is \(x\in[a,b]\). If \(a=np.inf, b=np.inf\) then the function computes \(a=\)
np.amin(x)(), \(b=\)np.amax(x)(). The cdf is computed as:\[\begin{split}cdf(x, a, b) = \begin{cases} 0 \quad , x<a \\ \frac{x-a}{b-a} \quad , a\leq x\leq b \\ 1 \quad , x>b \end{cases}\end{split}\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.uniform.nb_uniform_pdf(x: ndarray, a: float, b: float) ndarray#
Normalised uniform probability density function, w/ args: a, b. Its range of support is \(x\in[a,b]\). If \(a=np.inf, b=np.inf\) then the function computes \(a=\)
np.amin(x)(), \(b=\)np.amax(x)(). The pdf is computed as:\[\begin{split}pdf(x, a, b) = \begin{cases} \frac{1}{b-a} \quad , a\leq x\leq b \\ 0 \quad , \text{otherwise} \end{cases}\end{split}\]As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.uniform.nb_uniform_scaled_cdf(x: ndarray, a: float, b: float, area: float) ndarray#
Uniform cdf scaled by the area, used for extended binned fits As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- pygama.math.functions.uniform.nb_uniform_scaled_pdf(x: ndarray, a: float, b: float, area: float) ndarray#
Scaled uniform probability distribution, w/ args: a, b. As a Numba JIT function, it runs slightly faster than ‘out of the box’ functions.
- class pygama.math.functions.uniform.uniform_gen(momtype=1, a=None, b=None, xtol=1e-14, badvalue=None, name=None, longname=None, shapes=None, extradoc=None, seed=None)#
Bases:
pygama_continuous