Owl_stats
Statistics: random number generators, PDF and CDF functions, and hypothesis tests. The module also includes some basic statistical functions such as mean, variance, skew, and etc.
sem x
calculates the standard error of x
, also referred to as standard error of the mean.
absdev x
calculates the average absolute deviation of x
.
skew x
calculates the skewness (the third standardized moment) of x
.
kurtosis x
calculates the Pearson's kurtosis of x
, i.e. the fourth standardized moment of x
.
central_moment n x
calculates the n
th central moment of x
.
cov x0 x1
calculates the covariance of x0
and x1
, the mean of x0
and x1
can be specified by m0
and m1
respectively.
concordant x y
calculates the number of concordant pairs in the given arrays x
and y
. A pair of observations \((x_i, y_i)\) and \((x_j, y_j)\) is concordant if both elements in one pair are either greater than or less than both elements in the other pair, i.e., \((x_i > x_j) \land (y_i > y_j)\) or \((x_i < x_j) \land (y_i < y_j)\).
The function returns the count of such pairs.
discordant x y
calculates the number of discordant pairs in the given arrays x
and y
. A pair of observations \((x_i, y_i)\) and \((x_j, y_j)\) is discordant if one element in one pair is greater than the corresponding element in the other pair, and the other element is smaller, i.e., \((x_i > x_j) \land (y_i < y_j)\) or \((x_i < x_j) \land (y_i > y_j)\).
The function returns the count of such pairs.
corrcoef x y
calculates the Pearson correlation of x
and y
. Namely, corrcoef x y = cov(x, y) / (sigma_x * sigma_y)
.
kendall_tau x y
calculates the Kendall Tau correlation between x
and y
.
spearman_rho x y
calculates the Spearman Rho correlation between x
and y
.
autocorrelation ~lag x
calculates the autocorrelation of x
with the given lag
.
percentile x p
returns the p
percentile of the data x
. p
is between 0. and 100. x
does not need to be sorted beforehand.
quantile x p
returns the p
quantile of the data x
. x
does not need to be sorted beforehand. When computing several quantiles on the same data, it is more efficient to "pre-apply" the function, as in ``let f = quantile x in List.map f 0.1 ; 0.5
``, in which case the data is sorted only once.
Raises Invalid_argument if p
is not between 0 and 1.
first_quartile x
returns the first quartile of x
, i.e. 25 percentiles.
third_quartile x
returns the third quartile of x
, i.e. 75 percentiles.
interquartile x
returns the interquartile range of x
which is defined as the first quartile subtracted from the third quartile.
minmax_i x
returns the indices of both minimum and maximum in x
.
sort x
sorts the elements in the x
in increasing order if inc = true
, otherwise in decreasing order if inc=false
. By default, inc
is true
. Note a copy is returned, the original data is not modified.
argsort x
sorts the elements in x
and returns the indices mapping of the elements in the current array to their original position in x
.
The sorting is in increasing order if inc = true
, otherwise in decreasing order if inc=false
. By default, inc
is true
.
Computes sample's ranks.
The ranking order is from the smallest one to the largest. For example rank [|54.; 74.; 55.; 86.; 56.|]
returns [|1.; 4.; 2.; 5.; 3.|]
. Note that the ranking starts with one!
ties_strategy
controls which ranks are assigned to equal values:
Average
the mean of ranks should be assigned to each value. Default.Min
the minimum of ranks is assigned to each value.Max
the maximum of ranks is assigned to each value.type histogram = Owl_base_stats.histogram
Type for computed histograms, with optional weighted counts and normalized counts.
val histogram :
[ `Bins of float array | `N of int ] ->
?weights:float array ->
float array ->
histogram
histogram bins x
creates a histogram from values in x
. If bins matches `N n
it will construct n
equally spaced bins from the minimum to the maximum in x
. If bins matches `Bins b
, b
is taken as the sorted array of boundaries of adjacent bin intervals. Bin boundaries are taken as left-inclusive, right-exclusive, except for the last bin which is also right-inclusive. Values outside the bins are dropped silently.
histogram bins ~weights x
creates a weighted histogram with the given weights
which must match x
in length. The bare counts are also provided.
Returns a histogram including the n+1
bin boundaries, n
counts and weighted counts if applicable, but without normalisation.
val histogram_sorted :
[ `Bins of float array | `N of int ] ->
?weights:float array ->
float array ->
histogram
histogram_sorted bins x
is like histogram
but assumes that x
is sorted already. This increases efficiency if there are less bins than data. Undefined results if x
is not in fact sorted.
normalize hist
calculates a probability mass function using hist.weighted_counts
if present, otherwise using hist.counts
. The result is stored in the normalised_counts
field and sums to one.
normalize_density hist
calculates a probability density function using hist.weighted_counts
if present, otherwise using hist.counts
. The result is normalized as a density that is piecewise constant over the bin intervals. That is, the sum over density times corresponding bin width is one. If bins are infinitely wide, their density is 0 and the sum over width times density of all finite bins is the total weight in the finite bins. The result is stored in the density
field.
val pp_hist : Stdlib.Format.formatter -> histogram -> unit
Pretty-print summary information on a histogram record
ecdf x
returns (x',f)
which are the empirical cumulative distribution function f
of x
at points x'
. x'
is just x
sorted in increasing order with duplicates removed. The function does not support nan
values in the array x
.
z_score x
calculates the z score of a given array x
.
normalise_pdf arr
takes an array of floats arr
representing a probability density function (PDF) and returns a new array where the values are normalized so that the sum of the array equals 1.
tukey_fences ?k x
returns a tuple of the lower and upper boundaries for values that are not outliers. k
defaults to the standard coefficient of 1.5
. For first and third quartiles Q1
and `Q3`, the range is computed as follows:
(Q1 - k*(Q3-Q1), Q3 + k*(Q3-Q1))
val gaussian_kde :
?bandwidth:[ `Silverman | `Scott ] ->
?n_points:int ->
float array ->
float array * float array
gaussian_kde x
is a Gaussian kernel density estimator. The estimation of the pdf runs in `O(sample_size * n_points)`, and returns an array tuple (a, b)
where a
is a uniformly spaced points from the sample range at which the density function was estimated, and b
is the estimates at these points.
Bandwidth selection rules is as follows: * Silverman: use `rule-of-thumb` for choosing the bandwidth. It defaults to 0.9 * min(SD, IQR / 1.34) * n^-0.2
. * Scott: same as Silverman, but with a factor, equal to 1.06
.
The default bandwidth value is Scott
.
metropolis_hastings target_density initial_state num_samples
performs the Metropolis-Hastings algorithm to generate samples from a target distribution.
The Metropolis-Hastings algorithm is a Markov Chain Monte Carlo (MCMC) method that generates a sequence of samples from a probability distribution for which direct sampling is difficult. The algorithm uses a proposal distribution to explore the state space, accepting or rejecting proposed moves based on the ratio of the target densities.
gibbs_sampling conditional_sampler initial_state num_samples
performs Gibbs sampling to generate samples from a multivariate distribution.
Gibbs sampling is a Markov Chain Monte Carlo (MCMC) algorithm used to generate samples from a joint distribution when direct sampling is difficult. It works by iteratively sampling each variable from its conditional distribution, given the current values of all other variables.
Record type contains the result of a hypothesis test.
val pp_hypothesis : Stdlib.Format.formatter -> hypothesis -> unit
Pretty printer of hypothesis type
val z_test :
mu:float ->
sigma:float ->
?alpha:float ->
?side:tail ->
float array ->
hypothesis
z_test ~mu ~sigma ~alpha ~side x
returns a test decision for the null hypothesis that the data x
comes from a normal distribution with mean mu
and a standard deviation sigma
, using the z-test of alpha
significance level. The alternative hypothesis is that the mean is not mu
.
The result (h,p,z)
: h
is true
if the test rejects the null hypothesis at the alpha
significance level, and false
otherwise. p
is the p-value and z
is the z-score.
val t_test :
mu:float ->
?alpha:float ->
?side:tail ->
float array ->
hypothesis
t_test ~mu ~alpha ~side x
returns a test decision of one-sample t-test which is a parametric test of the location parameter when the population standard deviation is unknown. mu
is population mean, alpha
is the significance level.
val t_test_paired :
?alpha:float ->
?side:tail ->
float array ->
float array ->
hypothesis
t_test_paired ~alpha ~side x y
returns a test decision for the null hypothesis that the data in x – y
comes from a normal distribution with mean equal to zero and unknown variance, using the paired-sample t-test.
val t_test_unpaired :
?alpha:float ->
?side:tail ->
?equal_var:bool ->
float array ->
float array ->
hypothesis
t_test_unpaired ~alpha ~side ~equal_var x y
returns a test decision for the null hypothesis that the data in vectors x
and y
comes from independent random samples from normal distributions with equal means and equal but unknown variances, using the two-sample t-test. The alternative hypothesis is that the data in x
and y
comes from populations with unequal means.
equal_var
indicates whether two samples have the same variance. If the two variances are not the same, the test is referred to as Welche's t-test.
val ks_test : ?alpha:float -> float array -> (float -> float) -> hypothesis
ks_test ~alpha x f
returns a test decision for the null hypothesis that the data in vector x
comes from independent random samples of the distribution with CDF f. The alternative hypothesis is that the data in x
comes from a different distribution.
The result (h,p,d)
: h
is true
if the test rejects the null hypothesis at the alpha
significance level, and false
otherwise. p
is the p-value and d
is the Kolmogorov-Smirnov test statistic.
val ks2_test : ?alpha:float -> float array -> float array -> hypothesis
ks2_test ~alpha x y
returns a test decision for the null hypothesis that the data in vectors x
and y
come from independent random samples of the same distribution. The alternative hypothesis is that the data in x
and y
are sampled from different distributions.
The result (h,p,d)
: h
is true
if the test rejects the null hypothesis at the alpha
significance level, and false
otherwise. p
is the p-value and d
is the Kolmogorov-Smirnov test statistic.
val var_test :
?alpha:float ->
?side:tail ->
variance:float ->
float array ->
hypothesis
var_test ~alpha ~side ~variance x
returns a test decision for the null hypothesis that the data in x
comes from a normal distribution with input variance
, using the chi-square variance test. The alternative hypothesis is that x
comes from a normal distribution with a different variance.
val jb_test : ?alpha:float -> float array -> hypothesis
jb_test ~alpha x
returns a test decision for the null hypothesis that the data x
comes from a normal distribution with an unknown mean and variance, using the Jarque-Bera test.
val fisher_test :
?alpha:float ->
?side:tail ->
int ->
int ->
int ->
int ->
hypothesis
fisher_test ~alpha ~side a b c d
fisher's exact test for contingency table | a
, b
| | c
, d
|
The result (h,p,z)
: h
is true
if the test rejects the null hypothesis at the alpha
significance level, and false
otherwise. p
is the p-value and z
is prior odds ratio.
val runs_test :
?alpha:float ->
?side:tail ->
?v:float ->
float array ->
hypothesis
runs_test ~alpha ~v x
returns a test decision for the null hypothesis that the data x
comes in random order, against the alternative that they do not, by running Wald–Wolfowitz runs test. The test is based on the number of runs of consecutive values above or below the mean of x
. ~v
is the reference value, the default value is the median of x
.
val mannwhitneyu :
?alpha:float ->
?side:tail ->
float array ->
float array ->
hypothesis
mannwhitneyu ~alpha ~side x y
Computes the Mann-Whitney rank test on samples x and y. If length of each sample less than 10 and no ties, then using exact test (see paper Ying Kuen Cheung and Jerome H. Klotz (1997) The Mann Whitney Wilcoxon distribution using linked list Statistica Sinica 7 805-813), else usning asymptotic normal distribution.
val wilcoxon :
?alpha:float ->
?side:tail ->
float array ->
float array ->
hypothesis
wilcoxon ?alpha ?side x y
performs the Wilcoxon signed-rank test on the paired samples x
and y
.
The Wilcoxon signed-rank test is a non-parametric statistical hypothesis test used to compare two related samples, matched samples, or repeated measurements on a single sample to assess whether their population mean ranks differ.
The _rvs
functions generate random numbers according to the specified distribution. _pdf
are "density" functions that return the probability of the element specified by the arguments, while _cdf
functions are cumulative distribution functions that return the probability of all elements less than or equal to the chosen element, and _sf
functions are survival functions returning one minus the corresponding CDF function. `log` versions of functions return the result for the natural logarithm of a chosen element.
uniform_rvs ~a ~b
returns a random uniformly distributed integer between a
and b
, inclusive.
binomial_rvs p n
returns a random integer representing the number of successes in n
trials with probability of success p
on each trial.
binomial_pdf k ~p ~n
returns the binomially distributed probability of k
successes in n
trials with probability p
of success on each trial.
binomial_logpdf k ~p ~n
returns the log-binomially distributed probability of k
successes in n
trials with probability p
of success on each trial.
binomial_cdf k ~p ~n
returns the binomially distributed cumulative probability of less than or equal to k
successes in n
trials, with probability p
on each trial.
binomial_logcdf k ~p ~n
returns the log-binomially distributed cumulative probability of less than or equal to k
successes in n
trials, with probability p
on each trial.
binomial_sf k ~p ~n
is the binomial survival function, i.e. 1 - (binomial_cdf k ~p ~n)
.
binomial_loggf k ~p ~n
is the logbinomial survival function, i.e. 1 - (binomial_logcdf k ~p ~n)
.
hypergeometric_rvs ~good ~bad ~sample
returns a random hypergeometrically distributed integer representing the number of successes in a sample (without replacement) of size ~sample
from a population with ~good
successful elements and ~bad
unsuccessful elements.
hypergeometric_pdf k ~good ~bad ~sample
returns the hypergeometrically distributed probability of k
successes in a sample (without replacement) of ~sample
elements from a population containing ~good
successful elements and ~bad
unsuccessful ones.
hypergeometric_logpdf k ~good ~bad ~sample
returns a value equivalent to a log-transformed result from hypergeometric_pdf
.
multinomial_rvs n ~p
generates random numbers of multinomial distribution from n
trials. The probability mass function is as follows.
P(x) = \frac{n!}{{x_1}! \cdot\cdot\cdot {x_k}!} p_{1}^{x_1} \cdot\cdot\cdot p_{k}^{x_k}
p
is the probability mass of k
categories. The elements in p
should all be positive (result is undefined if there are negative values), but they don't need to sum to 1: the result is the same whether or not p
is normalized. For implementation details, refer to :cite:`davis1993computer`.
multinomial_rvs x ~p
return the probability of x
given the probability mass of a multinomial distribution.
multinomial_rvs x ~p
returns the logarithm probability of x
given the probability mass of a multinomial distribution.
categorical_rvs p
returns the value of a random variable which follows the categorical distribution. This is equavalent to only one trial from multinomial_rvs
function, so it is just a simple wrapping.
Continuous random variables
std_uniform_rvs ()
generates a random variate from the standard uniform distribution over the interval [0, 1\).
uniform_rvs ~a ~b
generates a random variate from the uniform distribution over the interval [a, b\).
uniform_pdf x ~a ~b
computes the probability density function (PDF) of the uniform distribution at the point x
over the interval [a, b\).
uniform_logpdf x ~a ~b
computes the natural logarithm of the probability density function (log-PDF) of the uniform distribution at the point x
over the interval [a, b\).
uniform_cdf x ~a ~b
computes the cumulative distribution function (CDF) of the uniform distribution at the point x
over the interval [a, b\).
uniform_logcdf x ~a ~b
computes the natural logarithm of the cumulative distribution function (log-CDF) of the uniform distribution at the point x
over the interval [a, b\).
uniform_ppf q ~a ~b
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the uniform distribution for a given probability q
over the interval [a, b\).
uniform_sf x ~a ~b
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the uniform distribution at the point x
over the interval [a, b\).
uniform_logsf x ~a ~b
computes the natural logarithm of the survival function (log-SF) of the uniform distribution at the point x
over the interval [a, b\).
uniform_isf q ~a ~b
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the uniform distribution for a given probability q
over the interval [a, b\).
exponential_rvs ~lambda
generates a random variate from the exponential distribution with rate parameter lambda
.
exponential_pdf x ~lambda
computes the probability density function (PDF) of the exponential distribution at the point x
with rate parameter lambda
.
exponential_logpdf x ~lambda
computes the natural logarithm of the probability density function (log-PDF) of the exponential distribution at the point x
with rate parameter lambda
.
exponential_cdf x ~lambda
computes the cumulative distribution function (CDF) of the exponential distribution at the point x
with rate parameter lambda
.
exponential_logcdf x ~lambda
computes the natural logarithm of the cumulative distribution function (log-CDF) of the exponential distribution at the point x
with rate parameter lambda
.
exponential_ppf q ~lambda
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the exponential distribution for a given probability q
with rate parameter lambda
.
exponential_sf x ~lambda
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the exponential distribution at the point x
with rate parameter lambda
.
exponential_logsf x ~lambda
computes the natural logarithm of the survival function (log-SF) of the exponential distribution at the point x
with rate parameter lambda
.
exponential_isf q ~lambda
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the exponential distribution for a given probability q
with rate parameter lambda
.
exponpow_pdf x ~a ~b
computes the probability density function (PDF) of the exponential power distribution at the point x
with shape parameters a
and b
.
exponpow_logpdf x ~a ~b
computes the natural logarithm of the probability density function (log-PDF) of the exponential power distribution at the point x
with shape parameters a
and b
.
exponpow_cdf x ~a ~b
computes the cumulative distribution function (CDF) of the exponential power distribution at the point x
with shape parameters a
and b
.
exponpow_logcdf x ~a ~b
computes the natural logarithm of the cumulative distribution function (log-CDF) of the exponential power distribution at the point x
with shape parameters a
and b
.
exponpow_sf x ~a ~b
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the exponential power distribution at the point x
with shape parameters a
and b
.
exponpow_logsf x ~a ~b
computes the natural logarithm of the survival function (log-SF) of the exponential power distribution at the point x
with shape parameters a
and b
.
gaussian_rvs ~mu ~sigma
generates a random variate from the Gaussian (normal) distribution with mean mu
and standard deviation sigma
.
gaussian_pdf x ~mu ~sigma
computes the probability density function (PDF) of the Gaussian (normal) distribution at the point x
with mean mu
and standard deviation sigma
.
gaussian_logpdf x ~mu ~sigma
computes the natural logarithm of the probability density function (log-PDF) of the Gaussian (normal) distribution at the point x
with mean mu
and standard deviation sigma
.
gaussian_cdf x ~mu ~sigma
computes the cumulative distribution function (CDF) of the Gaussian (normal) distribution at the point x
with mean mu
and standard deviation sigma
.
gaussian_logcdf x ~mu ~sigma
computes the natural logarithm of the cumulative distribution function (log-CDF) of the Gaussian (normal) distribution at the point x
with mean mu
and standard deviation sigma
.
gaussian_ppf q ~mu ~sigma
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the Gaussian (normal) distribution for a given probability q
with mean mu
and standard deviation sigma
.
gaussian_sf x ~mu ~sigma
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the Gaussian (normal) distribution at the point x
with mean mu
and standard deviation sigma
.
gaussian_logsf x ~mu ~sigma
computes the natural logarithm of the survival function (log-SF) of the Gaussian (normal) distribution at the point x
with mean mu
and standard deviation sigma
.
gaussian_isf q ~mu ~sigma
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the Gaussian (normal) distribution for a given probability q
with mean mu
and standard deviation sigma
.
gamma_rvs ~shape ~scale
generates a random variate from the gamma distribution with shape parameter shape
and scale parameter scale
.
gamma_pdf x ~shape ~scale
computes the probability density function (PDF) of the gamma distribution at the point x
with shape parameter shape
and scale parameter scale
.
gamma_logpdf x ~shape ~scale
computes the natural logarithm of the probability density function (log-PDF) of the gamma distribution at the point x
with shape parameter shape
and scale parameter scale
.
gamma_cdf x ~shape ~scale
computes the cumulative distribution function (CDF) of the gamma distribution at the point x
with shape parameter shape
and scale parameter scale
.
gamma_logcdf x ~shape ~scale
computes the natural logarithm of the cumulative distribution function (log-CDF) of the gamma distribution at the point x
with shape parameter shape
and scale parameter scale
.
gamma_ppf q ~shape ~scale
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the gamma distribution for a given probability q
with shape parameter shape
and scale parameter scale
.
gamma_sf x ~shape ~scale
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the gamma distribution at the point x
with shape parameter shape
and scale parameter scale
.
gamma_logsf x ~shape ~scale
computes the natural logarithm of the survival function (log-SF) of the gamma distribution at the point x
with shape parameter shape
and scale parameter scale
.
gamma_isf q ~shape ~scale
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the gamma distribution for a given probability q
with shape parameter shape
and scale parameter scale
.
beta_rvs ~a ~b
generates a random variate from the beta distribution with shape parameters a
and b
.
beta_pdf x ~a ~b
computes the probability density function (PDF) of the beta distribution at the point x
with shape parameters a
and b
.
beta_logpdf x ~a ~b
computes the natural logarithm of the probability density function (log-PDF) of the beta distribution at the point x
with shape parameters a
and b
.
beta_cdf x ~a ~b
computes the cumulative distribution function (CDF) of the beta distribution at the point x
with shape parameters a
and b
.
beta_logcdf x ~a ~b
computes the natural logarithm of the cumulative distribution function (log-CDF) of the beta distribution at the point x
with shape parameters a
and b
.
beta_ppf q ~a ~b
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the beta distribution for a given probability q
with shape parameters a
and b
.
beta_sf x ~a ~b
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the beta distribution at the point x
with shape parameters a
and b
.
beta_logsf x ~a ~b
computes the natural logarithm of the survival function (log-SF) of the beta distribution at the point x
with shape parameters a
and b
.
beta_isf q ~a ~b
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the beta distribution for a given probability q
with shape parameters a
and b
.
chi2_rvs ~df
generates a random variate from the chi-square distribution with degrees of freedom df
.
chi2_pdf x ~df
computes the probability density function (PDF) of the chi-square distribution at the point x
with degrees of freedom df
.
chi2_logpdf x ~df
computes the natural logarithm of the probability density function (log-PDF) of the chi-square distribution at the point x
with degrees of freedom df
.
chi2_cdf x ~df
computes the cumulative distribution function (CDF) of the chi-square distribution at the point x
with degrees of freedom df
.
chi2_logcdf x ~df
computes the natural logarithm of the cumulative distribution function (log-CDF) of the chi-square distribution at the point x
with degrees of freedom df
.
chi2_ppf q ~df
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the chi-square distribution for a given probability q
with degrees of freedom df
.
chi2_sf x ~df
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the chi-square distribution at the point x
with degrees of freedom df
.
chi2_logsf x ~df
computes the natural logarithm of the survival function (log-SF) of the chi-square distribution at the point x
with degrees of freedom df
.
chi2_isf q ~df
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the chi-square distribution for a given probability q
with degrees of freedom df
.
f_rvs ~dfnum ~dfden
generates a random variate from the F-distribution with numerator degrees of freedom dfnum
and denominator degrees of freedom dfden
.
The F-distribution is commonly used in the analysis of variance (ANOVA) and in the comparison of two variances.
f_pdf x ~dfnum ~dfden
computes the probability density function (PDF) of the F-distribution at the point x
with numerator degrees of freedom dfnum
and denominator degrees of freedom dfden
.
f_logpdf x ~dfnum ~dfden
computes the natural logarithm of the probability density function (log-PDF) of the F-distribution at the point x
with numerator degrees of freedom dfnum
and denominator degrees of freedom dfden
.
f_cdf x ~dfnum ~dfden
computes the cumulative distribution function (CDF) of the F-distribution at the point x
with numerator degrees of freedom dfnum
and denominator degrees of freedom dfden
.
f_logcdf x ~dfnum ~dfden
computes the natural logarithm of the cumulative distribution function (log-CDF) of the F-distribution at the point x
with numerator degrees of freedom dfnum
and denominator degrees of freedom dfden
.
f_ppf q ~dfnum ~dfden
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the F-distribution for a given probability q
with numerator degrees of freedom dfnum
and denominator degrees of freedom dfden
.
f_sf x ~dfnum ~dfden
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the F-distribution at the point x
with numerator degrees of freedom dfnum
and denominator degrees of freedom dfden
.
f_logsf x ~dfnum ~dfden
computes the natural logarithm of the survival function (log-SF) of the F-distribution at the point x
with numerator degrees of freedom dfnum
and denominator degrees of freedom dfden
.
f_isf q ~dfnum ~dfden
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the F-distribution for a given probability q
with numerator degrees of freedom dfnum
and denominator degrees of freedom dfden
.
cauchy_rvs ~loc ~scale
generates a random variate from the Cauchy distribution with location parameter loc
and scale parameter scale
.
The Cauchy distribution is a continuous probability distribution with heavy tails, often used in robust statistical methods.
cauchy_pdf x ~loc ~scale
computes the probability density function (PDF) of the Cauchy distribution at the point x
with location parameter loc
and scale parameter scale
.
cauchy_logpdf x ~loc ~scale
computes the natural logarithm of the probability density function (log-PDF) of the Cauchy distribution at the point x
with location parameter loc
and scale parameter scale
.
cauchy_cdf x ~loc ~scale
computes the cumulative distribution function (CDF) of the Cauchy distribution at the point x
with location parameter loc
and scale parameter scale
.
cauchy_logcdf x ~loc ~scale
computes the natural logarithm of the cumulative distribution function (log-CDF) of the Cauchy distribution at the point x
with location parameter loc
and scale parameter scale
.
cauchy_ppf q ~loc ~scale
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the Cauchy distribution for a given probability q
with location parameter loc
and scale parameter scale
.
cauchy_sf x ~loc ~scale
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the Cauchy distribution at the point x
with location parameter loc
and scale parameter scale
.
cauchy_logsf x ~loc ~scale
computes the natural logarithm of the survival function (log-SF) of the Cauchy distribution at the point x
with location parameter loc
and scale parameter scale
.
cauchy_isf q ~loc ~scale
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the Cauchy distribution for a given probability q
with location parameter loc
and scale parameter scale
.
t_rvs ~df ~loc ~scale
generates a random variate from the Student's t-distribution with df
degrees of freedom, location parameter loc
, and scale parameter scale
.
The Student's t-distribution is commonly used in statistics for small sample sizes or when the population standard deviation is unknown.
t_pdf x ~df ~loc ~scale
computes the probability density function (PDF) of the Student's t-distribution at the point x
with df
degrees of freedom, location parameter loc
, and scale parameter scale
.
t_logpdf x ~df ~loc ~scale
computes the natural logarithm of the probability density function (log-PDF) of the Student's t-distribution at the point x
with df
degrees of freedom, location parameter loc
, and scale parameter scale
.
t_cdf x ~df ~loc ~scale
computes the cumulative distribution function (CDF) of the Student's t-distribution at the point x
with df
degrees of freedom, location parameter loc
, and scale parameter scale
.
t_logcdf x ~df ~loc ~scale
computes the natural logarithm of the cumulative distribution function (log-CDF) of the Student's t-distribution at the point x
with df
degrees of freedom, location parameter loc
, and scale parameter scale
.
t_ppf q ~df ~loc ~scale
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the Student's t-distribution for a given probability q
with df
degrees of freedom, location parameter loc
, and scale parameter scale
.
t_sf x ~df ~loc ~scale
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the Student's t-distribution at the point x
with df
degrees of freedom, location parameter loc
, and scale parameter scale
.
t_logsf x ~df ~loc ~scale
computes the natural logarithm of the survival function (log-SF) of the Student's t-distribution at the point x
with df
degrees of freedom, location parameter loc
, and scale parameter scale
.
t_isf q ~df ~loc ~scale
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the Student's t-distribution for a given probability q
with df
degrees of freedom, location parameter loc
, and scale parameter scale
.
vonmises_rvs ~mu ~kappa
generates a random variate from the von Mises distribution with mean direction mu
and concentration parameter kappa
.
The von Mises distribution is often used as a circular analogue of the normal distribution for data measured in angles or on a circle.
vonmises_pdf x ~mu ~kappa
computes the probability density function (PDF) of the von Mises distribution at the point x
with mean direction mu
and concentration parameter kappa
.
vonmises_logpdf x ~mu ~kappa
computes the natural logarithm of the probability density function (log-PDF) of the von Mises distribution at the point x
with mean direction mu
and concentration parameter kappa
.
vonmises_cdf x ~mu ~kappa
computes the cumulative distribution function (CDF) of the von Mises distribution at the point x
with mean direction mu
and concentration parameter kappa
.
vonmises_logcdf x ~mu ~kappa
computes the natural logarithm of the cumulative distribution function (log-CDF) of the von Mises distribution at the point x
with mean direction mu
and concentration parameter kappa
.
vonmises_sf x ~mu ~kappa
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the von Mises distribution at the point x
with mean direction mu
and concentration parameter kappa
.
vonmises_logsf x ~mu ~kappa
computes the natural logarithm of the survival function (log-SF) of the von Mises distribution at the point x
with mean direction mu
and concentration parameter kappa
.
lomax_rvs ~shape ~scale
generates a random variate from the Lomax distribution, also known as the Pareto distribution of the second kind, with shape parameter shape
and scale parameter scale
.
The Lomax distribution is often used in survival analysis and heavy-tailed modeling.
lomax_pdf x ~shape ~scale
computes the probability density function (PDF) of the Lomax distribution at the point x
with shape parameter shape
and scale parameter scale
.
lomax_logpdf x ~shape ~scale
computes the natural logarithm of the probability density function (log-PDF) of the Lomax distribution at the point x
with shape parameter shape
and scale parameter scale
.
lomax_cdf x ~shape ~scale
computes the cumulative distribution function (CDF) of the Lomax distribution at the point x
with shape parameter shape
and scale parameter scale
.
lomax_logcdf x ~shape ~scale
computes the natural logarithm of the cumulative distribution function (log-CDF) of the Lomax distribution at the point x
with shape parameter shape
and scale parameter scale
.
lomax_ppf q ~shape ~scale
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the Lomax distribution for a given probability q
with shape parameter shape
and scale parameter scale
.
lomax_sf x ~shape ~scale
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the Lomax distribution at the point x
with shape parameter shape
and scale parameter scale
.
lomax_logsf x ~shape ~scale
computes the natural logarithm of the survival function (log-SF) of the Lomax distribution at the point x
with shape parameter shape
and scale parameter scale
.
lomax_isf q ~shape ~scale
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the Lomax distribution for a given probability q
with shape parameter shape
and scale parameter scale
.
weibull_rvs ~shape ~scale
generates a random variate from the Weibull distribution with shape parameter shape
and scale parameter scale
.
The Weibull distribution is commonly used in reliability analysis and modeling life data.
weibull_pdf x ~shape ~scale
computes the probability density function (PDF) of the Weibull distribution at the point x
with shape parameter shape
and scale parameter scale
.
weibull_logpdf x ~shape ~scale
computes the natural logarithm of the probability density function (log-PDF) of the Weibull distribution at the point x
with shape parameter shape
and scale parameter scale
.
weibull_cdf x ~shape ~scale
computes the cumulative distribution function (CDF) of the Weibull distribution at the point x
with shape parameter shape
and scale parameter scale
.
weibull_logcdf x ~shape ~scale
computes the natural logarithm of the cumulative distribution function (log-CDF) of the Weibull distribution at the point x
with shape parameter shape
and scale parameter scale
.
weibull_ppf q ~shape ~scale
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the Weibull distribution for a given probability q
with shape parameter shape
and scale parameter scale
.
weibull_sf x ~shape ~scale
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the Weibull distribution at the point x
with shape parameter shape
and scale parameter scale
.
weibull_logsf x ~shape ~scale
computes the natural logarithm of the survival function (log-SF) of the Weibull distribution at the point x
with shape parameter shape
and scale parameter scale
.
weibull_isf q ~shape ~scale
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the Weibull distribution for a given probability q
with shape parameter shape
and scale parameter scale
.
laplace_rvs ~loc ~scale
generates a random variate from the Laplace distribution, also known as the double exponential distribution, with location parameter loc
and scale parameter scale
.
The Laplace distribution is often used in statistical methods for detecting outliers.
laplace_pdf x ~loc ~scale
computes the probability density function (PDF) of the Laplace distribution at the point x
with location parameter loc
and scale parameter scale
.
laplace_logpdf x ~loc ~scale
computes the natural logarithm of the probability density function (log-PDF) of the Laplace distribution at the point x
with location parameter loc
and scale parameter scale
.
laplace_cdf x ~loc ~scale
computes the cumulative distribution function (CDF) of the Laplace distribution at the point x
with location parameter loc
and scale parameter scale
.
laplace_logcdf x ~loc ~scale
computes the natural logarithm of the cumulative distribution function (log-CDF) of the Laplace distribution at the point x
with location parameter loc
and scale parameter scale
.
laplace_ppf q ~loc ~scale
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the Laplace distribution for a given probability q
with location parameter loc
and scale parameter scale
.
laplace_sf x ~loc ~scale
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the Laplace distribution at the point x
with location parameter loc
and scale parameter scale
.
laplace_logsf x ~loc ~scale
computes the natural logarithm of the survival function (log-SF) of the Laplace distribution at the point x
with location parameter loc
and scale parameter scale
.
laplace_isf q ~loc ~scale
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the Laplace distribution for a given probability q
with location parameter loc
and scale parameter scale
.
gumbel1_rvs ~a ~b
generates a random variate from the Gumbel Type I distribution with location parameter a
and scale parameter b
.
The Gumbel distribution is commonly used to model the distribution of the maximum (or the minimum) of a number of samples of various distributions.
gumbel1_pdf x ~a ~b
computes the probability density function (PDF) of the Gumbel Type I distribution at the point x
with location parameter a
and scale parameter b
.
gumbel1_logpdf x ~a ~b
computes the natural logarithm of the probability density function (log-PDF) of the Gumbel Type I distribution at the point x
with location parameter a
and scale parameter b
.
gumbel1_cdf x ~a ~b
computes the cumulative distribution function (CDF) of the Gumbel Type I distribution at the point x
with location parameter a
and scale parameter b
.
gumbel1_logcdf x ~a ~b
computes the natural logarithm of the cumulative distribution function (log-CDF) of the Gumbel Type I distribution at the point x
with location parameter a
and scale parameter b
.
gumbel1_ppf q ~a ~b
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the Gumbel Type I distribution for a given probability q
with location parameter a
and scale parameter b
.
gumbel1_sf x ~a ~b
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the Gumbel Type I distribution at the point x
with location parameter a
and scale parameter b
.
gumbel1_logsf x ~a ~b
computes the natural logarithm of the survival function (log-SF) of the Gumbel Type I distribution at the point x
with location parameter a
and scale parameter b
.
gumbel1_isf q ~a ~b
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the Gumbel Type I distribution for a given probability q
with location parameter a
and scale parameter b
.
gumbel2_rvs ~a ~b
generates a random variate from the Gumbel Type II distribution with location parameter a
and scale parameter b
.
The Gumbel Type II distribution is used for modeling extreme values in various fields.
gumbel2_pdf x ~a ~b
computes the probability density function (PDF) of the Gumbel Type II distribution at the point x
with location parameter a
and scale parameter b
.
gumbel2_logpdf x ~a ~b
computes the natural logarithm of the probability density function (log-PDF) of the Gumbel Type II distribution at the point x
with location parameter a
and scale parameter b
.
gumbel2_cdf x ~a ~b
computes the cumulative distribution function (CDF) of the Gumbel Type II distribution at the point x
with location parameter a
and scale parameter b
.
gumbel2_logcdf x ~a ~b
computes the natural logarithm of the cumulative distribution function (log-CDF) of the Gumbel Type II distribution at the point x
with location parameter a
and scale parameter b
.
gumbel2_ppf q ~a ~b
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the Gumbel Type II distribution for a given probability q
with location parameter a
and scale parameter b
.
gumbel2_sf x ~a ~b
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the Gumbel Type II distribution at the point x
with location parameter a
and scale parameter b
.
gumbel2_logsf x ~a ~b
computes the natural logarithm of the survival function (log-SF) of the Gumbel Type II distribution at the point x
with location parameter a
and scale parameter b
.
gumbel2_isf q ~a ~b
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the Gumbel Type II distribution for a given probability q
with location parameter a
and scale parameter b
.
logistic_rvs ~loc ~scale
generates a random variate from the logistic distribution with location parameter loc
and scale parameter scale
.
The logistic distribution is often used in logistic regression and as a growth model.
logistic_pdf x ~loc ~scale
computes the probability density function (PDF) of the logistic distribution at the point x
with location parameter loc
and scale parameter scale
.
logistic_logpdf x ~loc ~scale
computes the natural logarithm of the probability density function (log-PDF) of the logistic distribution at the point x
with location parameter loc
and scale parameter scale
.
logistic_cdf x ~loc ~scale
computes the cumulative distribution function (CDF) of the logistic distribution at the point x
with location parameter loc
and scale parameter scale
.
logistic_logcdf x ~loc ~scale
computes the natural logarithm of the cumulative distribution function (log-CDF) of the logistic distribution at the point x
with location parameter loc
and scale parameter scale
.
logistic_ppf q ~loc ~scale
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the logistic distribution for a given probability q
with location parameter loc
and scale parameter scale
.
logistic_sf x ~loc ~scale
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the logistic distribution at the point x
with location parameter loc
and scale parameter scale
.
logistic_logsf x ~loc ~scale
computes the natural logarithm of the survival function (log-SF) of the logistic distribution at the point x
with location parameter loc
and scale parameter scale
.
logistic_isf q ~loc ~scale
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the logistic distribution for a given probability q
with location parameter loc
and scale parameter scale
.
lognormal_rvs ~mu ~sigma
generates a random variate from the log-normal distribution with parameters mu
(mean of the underlying normal distribution) and sigma
(standard deviation of the underlying normal distribution).
The log-normal distribution is commonly used to model positive-valued data with a distribution that is skewed to the right.
lognormal_pdf x ~mu ~sigma
computes the probability density function (PDF) of the log-normal distribution at the point x
with parameters mu
(mean of the underlying normal distribution) and sigma
(standard deviation of the underlying normal distribution).
lognormal_logpdf x ~mu ~sigma
computes the natural logarithm of the probability density function (log-PDF) of the log-normal distribution at the point x
with parameters mu
and sigma
.
lognormal_cdf x ~mu ~sigma
computes the cumulative distribution function (CDF) of the log-normal distribution at the point x
with parameters mu
and sigma
.
lognormal_logcdf x ~mu ~sigma
computes the natural logarithm of the cumulative distribution function (log-CDF) of the log-normal distribution at the point x
with parameters mu
and sigma
.
lognormal_ppf q ~mu ~sigma
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the log-normal distribution for a given probability q
with parameters mu
and sigma
.
lognormal_sf x ~mu ~sigma
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the log-normal distribution at the point x
with parameters mu
and sigma
.
lognormal_logsf x ~mu ~sigma
computes the natural logarithm of the survival function (log-SF) of the log-normal distribution at the point x
with parameters mu
and sigma
.
lognormal_isf q ~mu ~sigma
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the log-normal distribution for a given probability q
with parameters mu
and sigma
.
rayleigh_rvs ~sigma
generates a random variate from the Rayleigh distribution with scale parameter sigma
.
The Rayleigh distribution is commonly used to model the magnitude of a vector in two-dimensional space with Gaussian-distributed components.
rayleigh_pdf x ~sigma
computes the probability density function (PDF) of the Rayleigh distribution at the point x
with scale parameter sigma
.
rayleigh_logpdf x ~sigma
computes the natural logarithm of the probability density function (log-PDF) of the Rayleigh distribution at the point x
with scale parameter sigma
.
rayleigh_cdf x ~sigma
computes the cumulative distribution function (CDF) of the Rayleigh distribution at the point x
with scale parameter sigma
.
rayleigh_logcdf x ~sigma
computes the natural logarithm of the cumulative distribution function (log-CDF) of the Rayleigh distribution at the point x
with scale parameter sigma
.
rayleigh_ppf q ~sigma
computes the percent-point function (PPF), also known as the quantile function or inverse CDF, of the Rayleigh distribution for a given probability q
with scale parameter sigma
.
rayleigh_sf x ~sigma
computes the survival function (SF), which is one minus the cumulative distribution function (1 - CDF), of the Rayleigh distribution at the point x
with scale parameter sigma
.
rayleigh_logsf x ~sigma
computes the natural logarithm of the survival function (log-SF) of the Rayleigh distribution at the point x
with scale parameter sigma
.
rayleigh_isf q ~sigma
computes the inverse survival function (ISF), which is the inverse of the survival function (SF), of the Rayleigh distribution for a given probability q
with scale parameter sigma
.
dirichlet_rvs ~alpha
returns random variables of K-1
order Dirichlet distribution, follows the following probability dense function.
f(x_1,...,x_K; \alpha_1,...,\alpha_K) = \frac{1}{\mathbf{B(\alpha)}} \prod_{i=1}^K x_i^{\alpha_i - 1}
The normalising constant is the multivariate Beta function, which can be expressed in terms of the gamma function:
\mathbf{B(\alpha)} = \frac{\prod_{i=1}^K \Gamma(\alpha_i)}{\Gamma(\sum_{i=1}^K \alpha_i)}
Note that x
is a standard K-1 simplex, i.e. \sum_i^K x_i = 1
and x_i \ge 0, \forall x_i \in [1,K]
.
dirichlet_pdf x ~alpha
computes the probability density function (PDF) of the Dirichlet distribution for the input vector x
with concentration parameters alpha
.
The Dirichlet distribution is a multivariate generalization of the Beta distribution and is commonly used as a prior distribution in Bayesian statistics, particularly in the context of categorical data and multinomial distributions.