Introduction to distrib
=======================

Package distrib contains a set of functions for making probability computations on both discrete and continuous univariate models.

What follows is a short reminder of basic probabilistic related definitions.

Let f(x) be the DENSITY FUNCTION of an absolute continuous random variable X. The DISTRIBUTION FUNCTION is defined as

                            x
                           /
                           [
                    F(x) = I     f(u) du
                           ]
                           /
                            minf

which equals the probability PR(X <= X).

The MEAN value is a localization parameter and is defined as

                          inf
                         /
                         [
                E[X]  =  I   x f(x) dx
                         ]
                         /
                          minf

The VARIANCE is a measure of variation,

                      inf
                     /
                     [                    2
              V[X] = I     f(x) (x - E[X])  dx
                     ]
                     /
                      minf

which is a positive real number. The square root of the variance is the STANDARD DEVIATION, D[X]=sqrt(V[X]), and it is another measure of variation.

The SKEWNESS COEFFICIENT is a measure of non-symmetry,

                      inf
                     /
                 1   [                    3
       SK[X] = ----- I     f(x) (x - E[X])  dx
                   3 ]
               D[X]  /
                      minf

And the KURTOSIS COEFFICIENT measures the peakedness of the distribution,

                      inf
                     /
                 1   [                    4
       KU[X] = ----- I     f(x) (x - E[X])  dx - 3
                   4 ]
               D[X]  /
                      minf

If X is gaussian, KU[X]=0. In fact, both skewness and kurtosis are shape parameters used to measure the non--gaussianity of a distribution.

If the random variable X is discrete, the density, or PROBABILITY, function f(x) takes positive values within certain countable set of numbers x_i, and zero elsewhere. In this case, the distribution function is

                            ====
                            \
                     F(x) =  >    f(x )
                            /        i
                            ====
                           x <= x
                            i

The mean, variance, standard deviation, skewness coefficient and kurtosis coefficient take the form

                            ====
                            \
                     E[X] =  >  x  f(x ) ,
                            /    i    i
                            ====
                             x
                              i

                     ====
                     \                     2
             V[X] =   >    f(x ) (x - E[X])  ,
                     /        i    i
                     ====
                      x
                       i

                    D[X] = sqrt(V[X]),

                          ====
                   1      \                     3
       SK[X] =  -------    >    f(x ) (x - E[X])
                D[X]^3    /        i    i
                          ====
                           x
                            i

and

                          ====
                   1      \                     4
       KU[X] =  -------    >    f(x ) (x - E[X])   - 3 ,
                D[X]^4    /        i    i
                          ====
                           x
                            i

respectively.

Package distrib includes functions for simulating random variates. Some of these functions make use of optional variables indicating the algorithm to be used. The general inverse method (based on the fact that if U is an uniform random number in (0,1), then F^(-1)(U) is a random variate with distribution F) is implemented in most cases; this is a suboptimal method in terms of timing, but useful for comparing with other algorithms. In this example, the perandom_formance of algorithms ahrens_cheng and inverse for simulating chi-square variates are compared by means of their histograms:

     (%i1) load(distrib)$
     (%i2) load(descriptive)$
     (%i3) showtime: true$
     Evaluation took 0.00 seconds (0.00 elapsed) using 32 bytes.
     (%i4) random_chi2_algorithm: 'ahrens_cheng$
                                            histogram(random_chi2(10,500))$
     Evaluation took 0.00 seconds (0.00 elapsed) using 40 bytes.
     Evaluation took 0.69 seconds (0.71 elapsed) using 5.694 MB.
     (%i6) random_chi2_algorithm: 'inverse$ histogram(random_chi2(10,500))$
     Evaluation took 0.00 seconds (0.00 elapsed) using 32 bytes.
     Evaluation took 10.15 seconds (10.17 elapsed) using 322.098 MB.

In order to make visual comparisons among algorithms for a discrete variate, function barsplot of the descriptive package should be used.

Note that some work remains to be done, since these simulating functions are not yet checked by more rigurous goodness of fit tests.

Please, consult an introductory manual on probability and statistics for more information about all this mathematical stuff.

There is a naming convention in package distrib. Every function name has two parts, the first one makes reference to the function or parameter we want to calculate,

     Functions:
        Density function            (pdf_*)
        Distribution function       (cdf_*)
        Quantile                    (quantile_*)
        Mean                        (mean_*)
        Variance                    (var_*)
        Standard deviation          (std_*)
        Skewness coefficient        (skewness_*)
        Kurtosis coefficient        (kurtosis_*)
        Random variate              (random_*)

The second part is an explicit reference to the probabilistic model,

     Continuous distributions:
        Normal              (*normal)
        Student             (*student_t)
        Chi^2               (*chi2)
        F                   (*f)
        Exponential         (*exp)
        Lognormal           (*lognormal)
        Gamma               (*gamma)
        Beta                (*beta)
        Continuous uniform  (*continuous_uniform)
        Logistic            (*logistic)
        Pareto              (*pareto)
        Weibull             (*weibull)
        Rayleigh            (*rayleigh)
        Laplace             (*laplace)
        Cauchy              (*cauchy)
        Gumbel              (*gumbel)

     Discrete distributions:
        Binomial             (*binomial)
        Poisson              (*poisson)
        Bernoulli            (*bernoulli)
        Geometric            (*geometric)
        Discrete uniform     (*discrete_uniform)
        hypergeometric       (*hypergeometric)
        Negative binomial    (*negative_binomial)

For example, pdf_student_t(x,n) is the density function of the Student distribution with N degrees of freedom, std_pareto(a,b) is the standard deviation of the Pareto distribution with parameters A and B and kurtosis_poisson(m) is the kurtosis coefficient of the Poisson distribution with mean M.

In order to make use of package distrib you need first to load it by typing

     (%i1) load(distrib)$

For comments, bugs or suggestions, please contact the author at 'MARIO AT EDU DOT XUNTA DOT ES'.

Functions and Variables for continuous distributions
====================================================

 -- Function: pdf_normal (X,M,S)
     Returns the value at X of the density function of a Normal(m,s) random variable, with s>0. To make use of this function, write first load(distrib).

 -- Function: cdf_normal (X,M,S)
     Returns the value at X of the distribution function of a Normal(m,s) random variable, with s>0. This function is defined in terms of Maxima's built-in error function erf.

          (%i1) load (distrib)$
          (%i2) assume(s>0)$ cdf_normal(x,m,s);
                                       x - m
                                 erf(---------)
                                     sqrt(2) s    1
          (%o3)                  -------------- + -
                                       2          2

     See also erf.

 -- Function: quantile_normal (Q,M,S)
     Returns the Q-quantile of a Normal(m,s) random variable, with s>0; in other words, this is the inverse of cdf_normal. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_normal (M,S)
     Returns the mean of a Normal(m,s) random variable, with s>0, namely M. To make use of this function, write first load(distrib).

 -- Function: var_normal (M,S)
     Returns the variance of a Normal(m,s) random variable, with s>0, namely S^2. To make use of this function, write first load(distrib).

 -- Function: std_normal (M,S)
     Returns the standard deviation of a Normal(m,s) random variable, with s>0, namely S. To make use of this function, write first load(distrib).

 -- Function: skewness_normal (M,S)
     Returns the skewness coefficient of a Normal(m,s) random variable, with s>0, which is always equal to 0. To make use of this function, write first load(distrib).

 -- Function: kurtosis_normal (M,S)
     Returns the kurtosis coefficient of a Normal(m,s) random variable, with s>0, which is always equal to 0. To make use of this function, write first load(distrib).

 -- Option variable: random_normal_algorithm
     Default value: box_mueller

     This is the selected algorithm for simulating random normal variates. Implemented algorithms are box_mueller and inverse:

       * box_mueller, based on algorithm described in Knuth, D.E. (1981) SEMINUMERICAL ALGORITHMS. THE ART OF COMPUTER PROGRAMMING. Addison-Wesley.

       * inverse, based on the general inverse method.

     See also random_normal.

 -- Function: random_normal (M,S)
 -- Function: random_normal (M,S,N)
     Returns a Normal(m,s) random variate, with s>0. Calling random_normal with a third argument N, a random sample of size N will be simulated.

     There are two algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_normal_algorithm, which defaults to box_mueller.

     See also random_normal_algorithm. To make use of this function, write first load(distrib).

 -- Function: pdf_student_t (X,N)
     Returns the value at X of the density function of a Student random variable t(n), with n>0. To make use of this function, write first load(distrib).

 -- Function: cdf_student_t (X,N)
     Returns the value at X of the distribution function of a Student random variable t(n), with n>0. This function has no closed form and it is numerically computed if the global variable numer equals true,  otherwise it returns a nominal expression.

          (%i1) load (distrib)$
          (%i2) cdf_student_t(1/2, 7/3);
                                               1  7
          (%o2)                  cdf_student_t(-, -)
                                               2  3
          (%i3) %,numer;
          (%o3)                   .6698450596140417

 -- Function: quantile_student_t (Q,N)
     Returns the Q-quantile of a Student random variable t(n), with n>0; in other words, this is the inverse of cdf_student_t. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_student_t (N)
     Returns the mean of a Student random variable t(n), with n>0, which is always equal to 0. To make use of this function, write first load(distrib).

 -- Function: var_student_t (N)
     Returns the variance of a Student random variable t(n), with n>2.

          (%i1) load (distrib)$
          (%i2) assume(n>2)$  var_student_t(n);
                                          n
          (%o3)                         -----
                                        n - 2

 -- Function: std_student_t (N)
     Returns the standard deviation of a Student random variable t(n), with n>2. To make use of this function, write first load(distrib).

 -- Function: skewness_student_t (N)
     Returns the skewness coefficient of a Student random variable t(n), with n>3, which is always equal to 0. To make use of this function, write first load(distrib).

 -- Function: kurtosis_student_t (N)
     Returns the kurtosis coefficient of a Student random variable t(n), with n>4. To make use of this function, write first load(distrib).

 -- Option variable: random_student_t_algorithm
     Default value: ratio

     This is the selected algorithm for simulating random Student variates. Implemented algorithms are inverse and ratio:

       * inverse, based on the general inverse method.

       * ratio, based on the fact that if Z is a normal random variable N(0,1) and S^2 is chi square  random variable with N degrees of freedom, Chi^2(n), then

                                          Z
                                X = -------------
                                    /   2  \ 1/2
                                    |  S   |
                                    | ---  |
                                    \  n   /

          is a Student random variable with N degrees of freedom, t(n).

     See also random_student_t.

 -- Function: random_student_t (N)
 -- Function: random_student_t (N,M)
     Returns a Student random variate t(n), with n>0. Calling random_student_t with a second argument M, a random sample of size M will be simulated.

     There are two algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_student_t_algorithm, which defaults to ratio.

     See also random_student_t_algorithm. To make use of this function, write first load(distrib).

 -- Function: pdf_chi2 (X,N)
     Returns the value at X of the density function of a Chi-square random variable Chi^2(n), with n>0.

     The Chi^2(n) random variable is equivalent to the Gamma(n/2,2), therefore when Maxima has not enough information to get the result, a noun form based on the gamma density is returned.

          (%i1) load (distrib)$
          (%i2) pdf_chi2(x,n);
                                              n
          (%o2)                  pdf_gamma(x, -, 2)
                                              2
          (%i3) assume(x>0, n>0)$  pdf_chi2(x,n);
                                   n/2 - 1   - x/2
                                  x        %e
          (%o4)                   ----------------
                                    n/2       n
                                   2    gamma(-)
                                              2

 -- Function: cdf_chi2 (X,N)
     Returns the value at X of the distribution function of a Chi-square random variable Chi^2(n), with n>0.

     This function has no closed form and it is numerically computed if the global variable numer equals true,  otherwise it returns a nominal expression based on the gamma distribution, since the Chi^2(n) random variable is equivalent to the Gamma(n/2,2).

          (%i1) load (distrib)$
          (%i2) cdf_chi2(3,4);
          (%o2)                  cdf_gamma(3, 2, 2)
          (%i3) cdf_chi2(3,4),numer;
          (%o3)                   .4421745996289249

 -- Function: quantile_chi2 (Q,N)
     Returns the Q-quantile of a Chi-square random variable Chi^2(n), with n>0; in other words, this is the inverse of cdf_chi2. Argument Q must be an element of [0,1].

     This function has no closed form and it is numerically computed if the global variable numer equals true,  otherwise it returns a nominal expression based on the gamma quantile function, since the Chi^2(n) random variable is equivalent to the Gamma(n/2,2).

          (%i1) load (distrib)$
          (%i2) quantile_chi2(0.99,9);
          (%o2)                   21.66599433346194
          (%i3) quantile_chi2(0.99,n);
                                                  n
          (%o3)              quantile_gamma(0.99, -, 2)
                                                  2

 -- Function: mean_chi2 (N)
     Returns the mean of a Chi-square random variable Chi^2(n), with n>0.

     The Chi^2(n) random variable is equivalent to the Gamma(n/2,2), therefore when Maxima has not enough information to get the result, a noun form based on the gamma mean is returned.

          (%i1) load (distrib)$
          (%i2) mean_chi2(n);
                                             n
          (%o2)                   mean_gamma(-, 2)
                                             2
          (%i3) assume(n>0)$ mean_chi2(n);
          (%o4)                           n

 -- Function: var_chi2 (N)
     Returns the variance of a Chi-square random variable Chi^2(n), with n>0.

     The Chi^2(n) random variable is equivalent to the Gamma(n/2,2), therefore when Maxima has not enough information to get the result, a noun form based on the gamma variance is returned.

          (%i1) load (distrib)$
          (%i2) var_chi2(n);
                                             n
          (%o2)                    var_gamma(-, 2)
                                             2
          (%i3) assume(n>0)$ var_chi2(n);
          (%o4)                          2 n

 -- Function: std_chi2 (N)
     Returns the standard deviation of a Chi-square random variable Chi^2(n), with n>0.

     The Chi^2(n) random variable is equivalent to the Gamma(n/2,2), therefore when Maxima has not enough information to get the result, a noun form based on the gamma standard deviation is returned.

          (%i1) load (distrib)$
          (%i2) std_chi2(n);
                                             n
          (%o2)                    std_gamma(-, 2)
                                             2
          (%i3) assume(n>0)$ std_chi2(n);
          (%o4)                    sqrt(2) sqrt(n)

 -- Function: skewness_chi2 (N)
     Returns the skewness coefficient of a Chi-square random variable Chi^2(n), with n>0.

     The Chi^2(n) random variable is equivalent to the Gamma(n/2,2), therefore when Maxima has not enough information to get the result, a noun form based on the gamma skewness coefficient is returned.

          (%i1) load (distrib)$
          (%i2) skewness_chi2(n);
                                               n
          (%o2)                 skewness_gamma(-, 2)
                                               2
          (%i3) assume(n>0)$ skewness_chi2(n);
                                      2 sqrt(2)
          (%o4)                       ---------
                                       sqrt(n)

 -- Function: kurtosis_chi2 (N)
     Returns the kurtosis coefficient of a Chi-square random variable Chi^2(n), with n>0.

     The Chi^2(n) random variable is equivalent to the Gamma(n/2,2), therefore when Maxima has not enough information to get the result, a noun form based on the gamma kurtosis coefficient is returned.

          (%i1) load (distrib)$
          (%i2) kurtosis_chi2(n);
                                               n
          (%o2)                 kurtosis_gamma(-, 2)
                                               2
          (%i3) assume(n>0)$ kurtosis_chi2(n);
                                         12
          (%o4)                          --
                                         n

 -- Option variable: random_chi2_algorithm
     Default value: ahrens_cheng

     This is the selected algorithm for simulating random Chi-square variates. Implemented algorithms are ahrens_cheng and inverse:

       * ahrens_cheng, based on the random simulation of gamma variates. See random_gamma_algorithm for details.

       * inverse, based on the general inverse method.

     See also random_chi2.

 -- Function: random_chi2 (N)
 -- Function: random_chi2 (N,M)
     Returns a Chi-square random variate Chi^2(n), with n>0. Calling random_chi2 with a second argument M, a random sample of size M will be simulated.

     There are two algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_chi2_algorithm, which defaults to ahrens_cheng.

     See also random_chi2_algorithm. To make use of this function, write first load(distrib).

 -- Function: pdf_f (X,M,N)
     Returns the value at X of the density function of a F random variable F(m,n), with m,n>0. To make use of this function, write first load(distrib).

 -- Function: cdf_f (X,M,N)
     Returns the value at X of the distribution function of a F random variable F(m,n), with m,n>0. This function has no closed form and it is numerically computed if the global variable numer equals true,  otherwise it returns a nominal expression.

          (%i1) load (distrib)$
          (%i2) cdf_f(2,3,9/4);
                                               9
          (%o2)                    cdf_f(2, 3, -)
                                               4
          (%i3) %,numer;
          (%o3)                   0.66756728179008

 -- Function: quantile_f (Q,M,N)
     Returns the Q-quantile of a F random variable F(m,n), with m,n>0; in other words, this is the inverse of cdf_f. Argument Q must be an element of [0,1].

     This function has no closed form and it is numerically computed if the global variable numer equals true,  otherwise it returns a nominal expression.

          (%i1) load (distrib)$
          (%i2) quantile_f(2/5,sqrt(3),5);
                                         2
          (%o2)               quantile_f(-, sqrt(3), 5)
                                         5
          (%i3) %,numer;
          (%o3)                   0.518947838573693

 -- Function: mean_f (M,N)
     Returns the mean of a F random variable F(m,n), with m>0, n>2. To make use of this function, write first load(distrib).

 -- Function: var_f (M,N)
     Returns the variance of a F random variable F(m,n), with m>0, n>4. To make use of this function, write first load(distrib).

 -- Function: std_f (M,N)
     Returns the standard deviation of a F random variable F(m,n), with m>0, n>4. To make use of this function, write first load(distrib).

 -- Function: skewness_f (M,N)
     Returns the skewness coefficient of a F random variable F(m,n), with m>0, n>6. To make use of this function, write first load(distrib).

 -- Function: kurtosis_f (M,N)
     Returns the kurtosis coefficient of a F random variable F(m,n), with m>0, n>8. To make use of this function, write first load(distrib).

 -- Option variable: random_f_algorithm
     Default value: inverse

     This is the selected algorithm for simulating random F variates. Implemented algorithms are ratio and inverse:

       * ratio, based on the fact that if X is a Chi^2(m) random variable and Y is a Chi^2(n) random variable, then

                                       n X
                                   F = ---
                                       m Y

          is a F random variable with M and N degrees of freedom, F(m,n).

       * inverse, based on the general inverse method.

     See also random_f.

 -- Function: random_f (M,N)
 -- Function: random_f (M,N,K)
     Returns a F random variate F(m,n), with m,n>0. Calling random_f with a third argument K, a random sample of size K will be simulated.

     There are two algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_f_algorithm, which defaults to inverse.

     See also random_f_algorithm. To make use of this function, write first load(distrib).

 -- Function: pdf_exp (X,M)
     Returns the value at X of the density function of an Exponential(m) random variable, with m>0.

     The Exponential(m) random variable is equivalent to the Weibull(1,1/m), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull density is returned.

          (%i1) load (distrib)$
          (%i2) pdf_exp(x,m);
                                                  1
          (%o2)                 pdf_weibull(x, 1, -)
                                                  m
          (%i3) assume(x>0,m>0)$  pdf_exp(x,m);
                                          - m x
          (%o4)                       m %e

 -- Function: cdf_exp (X,M)
     Returns the value at X of the distribution function of an Exponential(m) random variable, with m>0.

     The Exponential(m) random variable is equivalent to the Weibull(1,1/m), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull distribution is returned.

          (%i1) load (distrib)$
          (%i2) cdf_exp(x,m);
                                                  1
          (%o2)                 cdf_weibull(x, 1, -)
                                                  m
          (%i3) assume(x>0,m>0)$  cdf_exp(x,m);
                                           - m x
          (%o4)                      1 - %e

 -- Function: quantile_exp (Q,M)
     Returns the Q-quantile of an Exponential(m) random variable, with m>0; in other words, this is the inverse of cdf_exp. Argument Q must be an element of [0,1].

     The Exponential(m) random variable is equivalent to the Weibull(1,1/m), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull quantile is returned.

          (%i1) load (distrib)$
          (%i2) quantile_exp(0.56,5);
          (%o2)                   .1641961104139661
          (%i3) quantile_exp(0.56,m);
                                                      1
          (%o3)             quantile_weibull(0.56, 1, -)
                                                      m

 -- Function: mean_exp (M)
     Returns the mean of an Exponential(m) random variable, with m>0.

     The Exponential(m) random variable is equivalent to the Weibull(1,1/m), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull mean is returned.

          (%i1) load (distrib)$
          (%i2) mean_exp(m);
                                                 1
          (%o2)                  mean_weibull(1, -)
                                                 m
          (%i3) assume(m>0)$  mean_exp(m);
                                          1
          (%o4)                           -
                                          m

 -- Function: var_exp (M)
     Returns the variance of an Exponential(m) random variable, with m>0.

     The Exponential(m) random variable is equivalent to the Weibull(1,1/m), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull variance is returned.

          (%i1) load (distrib)$
          (%i2) var_exp(m);
                                                 1
          (%o2)                   var_weibull(1, -)
                                                 m
          (%i3) assume(m>0)$  var_exp(m);
                                         1
          (%o4)                          --
                                          2
                                         m

 -- Function: std_exp (M)
     Returns the standard deviation of an Exponential(m) random variable, with m>0.

     The Exponential(m) random variable is equivalent to the Weibull(1,1/m), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull standard deviation is returned.

          (%i1) load (distrib)$
          (%i2) std_exp(m);
                                                 1
          (%o2)                   std_weibull(1, -)
                                                 m
          (%i3) assume(m>0)$  std_exp(m);
                                          1
          (%o4)                           -
                                          m

 -- Function: skewness_exp (M)
     Returns the skewness coefficient of an Exponential(m) random variable, with m>0.

     The Exponential(m) random variable is equivalent to the Weibull(1,1/m), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull skewness coefficient is returned.

          (%i1) load (distrib)$
          (%i2) skewness_exp(m);
                                                   1
          (%o2)                skewness_weibull(1, -)
                                                   m
          (%i3) assume(m>0)$  skewness_exp(m);
          (%o4)                           2

 -- Function: kurtosis_exp (M)
     Returns the kurtosis coefficient of an Exponential(m) random variable, with m>0.

     The Exponential(m) random variable is equivalent to the Weibull(1,1/m), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull kurtosis coefficient is returned.

          (%i1) load (distrib)$
          (%i2) kurtosis_exp(m);
                                                   1
          (%o2)                kurtosis_weibull(1, -)
                                                   m
          (%i3) assume(m>0)$  kurtosis_exp(m);
          (%o4)                           6

 -- Option variable: random_exp_algorithm
     Default value: inverse

     This is the selected algorithm for simulating random exponential variates. Implemented algorithms are inverse, ahrens_cheng and ahrens_dieter

       * inverse, based on the general inverse method.

       * ahrens_cheng, based on the fact that the Exp(m) random variable is equivalent to the Gamma(1,1/m). See random_gamma_algorithm for details.

       * ahrens_dieter, based on algorithm described in Ahrens, J.H. and Dieter, U. (1972) COMPUTER METHODS FOR SAMPLING FROM THE EXPONENTIAL AND NORMAL DISTRIBUTIONS. Comm, ACM, 15, Oct.,  873-882.

     See also random_exp.

 -- Function: random_exp (M)
 -- Function: random_exp (M,K)
     Returns an Exponential(m) random variate, with m>0. Calling random_exp with a second argument K, a random sample of size K will be simulated.

     There are three algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_exp_algorithm, which defaults to inverse.

     See also random_exp_algorithm. To make use of this function, write first load(distrib).

 -- Function: pdf_lognormal (X,M,S)
     Returns the value at X of the density function of a Lognormal(m,s) random variable, with s>0. To make use of this function, write first load(distrib).

 -- Function: cdf_lognormal (X,M,S)
     Returns the value at X of the distribution function of a Lognormal(m,s) random variable, with s>0. This function is defined in terms of Maxima's built-in error function erf.

          (%i1) load (distrib)$
          (%i2) assume(x>0, s>0)$  cdf_lognormal(x,m,s);
                                     log(x) - m
                                 erf(----------)
                                     sqrt(2) s     1
          (%o3)                  --------------- + -
                                        2          2

     See also erf.

 -- Function: quantile_lognormal (Q,M,S)
     Returns the Q-quantile of a Lognormal(m,s) random variable, with s>0; in other words, this is the inverse of cdf_lognormal. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_lognormal (M,S)
     Returns the mean of a Lognormal(m,s) random variable, with s>0. To make use of this function, write first load(distrib).

 -- Function: var_lognormal (M,S)
     Returns the variance of a Lognormal(m,s) random variable, with s>0. To make use of this function, write first load(distrib).

 -- Function: std_lognormal (M,S)
     Returns the standard deviation of a Lognormal(m,s) random variable, with s>0. To make use of this function, write first load(distrib).

 -- Function: skewness_lognormal (M,S)
     Returns the skewness coefficient of a Lognormal(m,s) random variable, with s>0. To make use of this function, write first load(distrib).

 -- Function: kurtosis_lognormal (M,S)
     Returns the kurtosis coefficient of a Lognormal(m,s) random variable, with s>0. To make use of this function, write first load(distrib).

 -- Function: random_lognormal (M,S)
 -- Function: random_lognormal (M,S,N)
     Returns a Lognormal(m,s) random variate, with s>0. Calling random_lognormal with a third argument N, a random sample of size N will be simulated.

     Log-normal variates are simulated by means of random normal variates. There are two algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_normal_algorithm, which defaults to box_mueller.

     See also random_normal_algorithm. To make use of this function, write first load(distrib).

 -- Function: pdf_gamma (X,A,B)
     Returns the value at X of the density function of a Gamma(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: cdf_gamma (X,A,B)
     Returns the value at X of the distribution function of a Gamma(a,b) random variable, with a,b>0.

     This function has no closed form and it is numerically computed if the global variable numer equals true,  otherwise it returns a nominal expression.

          (%i1) load (distrib)$
          (%i2) cdf_gamma(3,5,21);
          (%o2)                  cdf_gamma(3, 5, 21)
          (%i3) %,numer;
          (%o3)                 4.402663157135039E-7

 -- Function: quantile_gamma (Q,A,B)
     Returns the Q-quantile of a Gamma(a,b) random variable, with a,b>0; in other words, this is the inverse of cdf_gamma. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_gamma (A,B)
     Returns the mean of a Gamma(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: var_gamma (A,B)
     Returns the variance of a Gamma(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: std_gamma (A,B)
     Returns the standard deviation of a Gamma(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: skewness_gamma (A,B)
     Returns the skewness coefficient of a Gamma(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: kurtosis_gamma (A,B)
     Returns the kurtosis coefficient of a Gamma(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Option variable: random_gamma_algorithm
     Default value: ahrens_cheng

     This is the selected algorithm for simulating random gamma variates. Implemented algorithms are ahrens_cheng and inverse

       * ahrens_cheng, this is a combinantion of two procedures, depending on the value of parameter A:

          For a>=1, Cheng, R.C.H. and Feast, G.M. (1979). SOME SIMPLE GAMMA VARIATE GENERATORS. Appl. Stat., 28, 3, 290-295.

          For 0<a<1, Ahrens, J.H. and Dieter, U. (1974). COMPUTER METHODS FOR SAMPLING FROM GAMMA, BETA, POISSON AND BINOMIAL CDF_TRIBUTIONS. Computing, 12, 223-246.

       * inverse, based on the general inverse method.

     See also random_gamma.

 -- Function: random_gamma (A,B)
 -- Function: random_gamma (A,B,N)
     Returns a Gamma(a,b) random variate, with a,b>0. Calling random_gamma with a third argument N, a random sample of size N will be simulated.

     There are two algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_gamma_algorithm, which defaults to ahrens_cheng.

     See also random_gamma_algorithm. To make use of this function, write first load(distrib).

 -- Function: pdf_beta (X,A,B)
     Returns the value at X of the density function of a Beta(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: cdf_beta (X,A,B)
     Returns the value at X of the distribution function of a Beta(a,b) random variable, with a,b>0.

     This function has no closed form and it is numerically computed if the global variable numer equals true,  otherwise it returns a nominal expression.

          (%i1) load (distrib)$
          (%i2) cdf_beta(1/3,15,2);
                                           1
          (%o2)                   cdf_beta(-, 15, 2)
                                           3
          (%i3) %,numer;
          (%o3)                 7.666089131388224E-7

 -- Function: quantile_beta (Q,A,B)
     Returns the Q-quantile of a Beta(a,b) random variable, with a,b>0; in other words, this is the inverse of cdf_beta. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_beta (A,B)
     Returns the mean of a Beta(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: var_beta (A,B)
     Returns the variance of a Beta(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: std_beta (A,B)
     Returns the standard deviation of a Beta(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: skewness_beta (A,B)
     Returns the skewness coefficient of a Beta(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: kurtosis_beta (A,B)
     Returns the kurtosis coefficient of a Beta(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Option variable: random_beta_algorithm
     Default value: cheng

     This is the selected algorithm for simulating random beta variates. Implemented algorithms are cheng, inverse and ratio

       * cheng, this is the algorithm defined in Cheng, R.C.H.  (1978). GENERATING BETA VARIATES WITH NONINTEGRAL SHAPE PARAMETERS. Communications of the ACM, 21:317-322

       * inverse, based on the general inverse method.

       * ratio, based on the fact that if X is a random variable Gamma(a,1) and Y is Gamma(b,1), then the ratio X/(X+Y) is distributed as Beta(a,b).

     See also random_beta.

 -- Function: random_beta (A,B)
 -- Function: random_beta (A,B,N)
     Returns a Beta(a,b) random variate, with a,b>0. Calling random_beta with a third argument N, a random sample of size N will be simulated.

     There are three algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_beta_algorithm, which defaults to cheng.

     See also random_beta_algorithm. To make use of this function, write first load(distrib).

 -- Function: pdf_continuous_uniform (X,A,B)
     Returns the value at X of the density function of a Continuous Uniform(a,b) random variable, with a<b. To make use of this function, write first load(distrib).

 -- Function: cdf_continuous_uniform (X,A,B)
     Returns the value at X of the distribution function of a Continuous Uniform(a,b) random variable, with a<b. To make use of this function, write first load(distrib).

 -- Function: quantile_continuous_uniform (Q,A,B)
     Returns the Q-quantile of a Continuous Uniform(a,b) random variable, with a<b; in other words, this is the inverse of cdf_continuous_uniform. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_continuous_uniform (A,B)
     Returns the mean of a Continuous Uniform(a,b) random variable, with a<b. To make use of this function, write first load(distrib).

 -- Function: var_continuous_uniform (A,B)
     Returns the variance of a Continuous Uniform(a,b) random variable, with a<b. To make use of this function, write first load(distrib).

 -- Function: std_continuous_uniform (A,B)
     Returns the standard deviation of a Continuous Uniform(a,b) random variable, with a<b. To make use of this function, write first load(distrib).

 -- Function: skewness_continuous_uniform (A,B)
     Returns the skewness coefficient of a Continuous Uniform(a,b) random variable, with a<b. To make use of this function, write first load(distrib).

 -- Function: kurtosis_continuous_uniform (A,B)
     Returns the kurtosis coefficient of a Continuous Uniform(a,b) random variable, with a<b. To make use of this function, write first load(distrib).

 -- Function: random_continuous_uniform (A,B)
 -- Function: random_continuous_uniform (A,B,N)
     Returns a Continuous Uniform(a,b) random variate, with a<b. Calling random_continuous_uniform with a third argument N, a random sample of size N will be simulated.

     This is a direct application of the random built-in Maxima function.

     See also random. To make use of this function, write first load(distrib).

 -- Function: pdf_logistic (X,A,B)
     Returns the value at X of the density function of a Logistic(a,b) random variable , with b>0. To make use of this function, write first load(distrib).

 -- Function: cdf_logistic (X,A,B)
     Returns the value at X of the distribution function of a Logistic(a,b) random variable , with b>0. To make use of this function, write first load(distrib).

 -- Function: quantile_logistic (Q,A,B)
     Returns the Q-quantile of a Logistic(a,b) random variable , with b>0; in other words, this is the inverse of cdf_logistic. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: meanlog (A,B)
     Returns the mean of a Logistic(a,b) random variable , with b>0. To make use of this function, write first load(distrib).

 -- Function: var_logistic (A,B)
     Returns the variance of a Logistic(a,b) random variable , with b>0. To make use of this function, write first load(distrib).

 -- Function: std_logistic (A,B)
     Returns the standard deviation of a Logistic(a,b) random variable , with b>0. To make use of this function, write first load(distrib).

 -- Function: skewness_logistic (A,B)
     Returns the skewness coefficient of a Logistic(a,b) random variable , with b>0. To make use of this function, write first load(distrib).

 -- Function: kurtosis_logistic (A,B)
     Returns the kurtosis coefficient of a Logistic(a,b) random variable , with b>0. To make use of this function, write first load(distrib).

 -- Function: random_logistic (A,B)
 -- Function: random_logistic (A,B,N)
     Returns a Logistic(a,b) random variate, with b>0. Calling random_logistic with a third argument N, a random sample of size N will be simulated.

     Only the inverse method is implemented. To make use of this function, write first load(distrib).

 -- Function: pdf_pareto (X,A,B)
     Returns the value at X of the density function of a Pareto(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: cdf_pareto (X,A,B)
     Returns the value at X of the distribution function of a Pareto(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: quantile_pareto (Q,A,B)
     Returns the Q-quantile of a Pareto(a,b) random variable, with a,b>0; in other words, this is the inverse of cdf_pareto. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_pareto (A,B)
     Returns the mean of a Pareto(a,b) random variable, with a>1,b>0. To make use of this function, write first load(distrib).

 -- Function: var_pareto (A,B)
     Returns the variance of a Pareto(a,b) random variable, with a>2,b>0. To make use of this function, write first load(distrib).

 -- Function: std_pareto (A,B)
     Returns the standard deviation of a Pareto(a,b) random variable, with a>2,b>0. To make use of this function, write first load(distrib).

 -- Function: skewness_pareto (A,B)
     Returns the skewness coefficient of a Pareto(a,b) random variable, with a>3,b>0. To make use of this function, write first load(distrib).

 -- Function: kurtosis_pareto (A,B)
     Returns the kurtosis coefficient of a Pareto(a,b) random variable, with a>4,b>0. To make use of this function, write first load(distrib).

 -- Function: random_pareto (A,B)
 -- Function: random_pareto (A,B,N)
     Returns a Pareto(a,b) random variate, with a>0,b>0. Calling random_pareto with a third argument N, a random sample of size N will be simulated.

     Only the inverse method is implemented. To make use of this function, write first load(distrib).

 -- Function: pdf_weibull (X,A,B)
     Returns the value at X of the density function of a Weibull(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: cdf_weibull (X,A,B)
     Returns the value at X of the distribution function of a Weibull(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: quantile_weibull (Q,A,B)
     Returns the Q-quantile of a Weibull(a,b) random variable, with a,b>0; in other words, this is the inverse of cdf_weibull. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_weibull (A,B)
     Returns the mean of a Weibull(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: var_weibull (A,B)
     Returns the variance of a Weibull(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: std_weibull (A,B)
     Returns the standard deviation of a Weibull(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: skewness_weibull (A,B)
     Returns the skewness coefficient of a Weibull(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: kurtosis_weibull (A,B)
     Returns the kurtosis coefficient of a Weibull(a,b) random variable, with a,b>0. To make use of this function, write first load(distrib).

 -- Function: random_weibull (A,B)
 -- Function: random_weibull (A,B,N)
     Returns a Weibull(a,b) random variate, with a,b>0. Calling random_weibull with a third argument N, a random sample of size N will be simulated.

     Only the inverse method is implemented. To make use of this function, write first load(distrib).

 -- Function: pdf_rayleigh (X,B)
     Returns the value at X of the density function of a Rayleigh(b) random variable, with b>0.

     The Rayleigh(b) random variable is equivalent to the Weibull(2,1/b), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull density is returned.

          (%i1) load (distrib)$
          (%i2) pdf_rayleigh(x,b);
                                                  1
          (%o2)                 pdf_weibull(x, 2, -)
                                                  b
          (%i3) assume(x>0,b>0)$ pdf_rayleigh(x,b);
                                              2  2
                                     2     - b  x
          (%o4)                   2 b  x %e

 -- Function: cdf_rayleigh (X,B)
     Returns the value at X of the distribution function of a Rayleigh(b) random variable, with b>0.

     The Rayleigh(b) random variable is equivalent to the Weibull(2,1/b), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull distribution is returned.

          (%i1) load (distrib)$
          (%i2) cdf_rayleigh(x,b);
                                                  1
          (%o2)                 cdf_weibull(x, 2, -)
                                                  b
          (%i3) assume(x>0,b>0)$ cdf_rayleigh(x,b);
                                             2  2
                                          - b  x
          (%o4)                     1 - %e

 -- Function: quantile_rayleigh (Q,B)
     Returns the Q-quantile of a Rayleigh(b) random variable, with b>0; in other words, this is the inverse of cdf_rayleigh. Argument Q must be an element of [0,1].

     The Rayleigh(b) random variable is equivalent to the Weibull(2,1/b), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull quantile is returned.

          (%i1) load (distrib)$
          (%i2) quantile_rayleigh(0.99,b);
                                                      1
          (%o2)             quantile_weibull(0.99, 2, -)
                                                      b
          (%i3) assume(x>0,b>0)$ quantile_rayleigh(0.99,b);
                                  2.145966026289347
          (%o4)                   -----------------
                                          b

 -- Function: mean_rayleigh (B)
     Returns the mean of a Rayleigh(b) random variable, with b>0.

     The Rayleigh(b) random variable is equivalent to the Weibull(2,1/b), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull mean is returned.

          (%i1) load (distrib)$
          (%i2) mean_rayleigh(b);
                                                 1
          (%o2)                  mean_weibull(2, -)
                                                 b
          (%i3) assume(b>0)$ mean_rayleigh(b);
                                      sqrt(%pi)
          (%o4)                       ---------
                                         2 b

 -- Function: var_rayleigh (B)
     Returns the variance of a Rayleigh(b) random variable, with b>0.

     The Rayleigh(b) random variable is equivalent to the Weibull(2,1/b), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull variance is returned.

          (%i1) load (distrib)$
          (%i2) var_rayleigh(b);
                                                 1
          (%o2)                   var_weibull(2, -)
                                                 b
          (%i3) assume(b>0)$ var_rayleigh(b);
                                           %pi
                                       1 - ---
                                            4
          (%o4)                        -------
                                          2
                                         b

 -- Function: std_rayleigh (B)
     Returns the standard deviation of a Rayleigh(b) random variable, with b>0.

     The Rayleigh(b) random variable is equivalent to the Weibull(2,1/b), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull standard deviation is returned.

          (%i1) load (distrib)$
          (%i2) std_rayleigh(b);
                                                 1
          (%o2)                   std_weibull(2, -)
                                                 b
          (%i3) assume(b>0)$ std_rayleigh(b);
                                             %pi
                                    sqrt(1 - ---)
                                              4
          (%o4)                     -------------
                                          b

 -- Function: skewness_rayleigh (B)
     Returns the skewness coefficient of a Rayleigh(b) random variable, with b>0.

     The Rayleigh(b) random variable is equivalent to the Weibull(2,1/b), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull skewness coefficient is returned.

          (%i1) load (distrib)$
          (%i2) skewness_rayleigh(b);
                                                   1
          (%o2)                skewness_weibull(2, -)
                                                   b
          (%i3) assume(b>0)$ skewness_rayleigh(b);
                                   3/2
                                %pi      3 sqrt(%pi)
                                ------ - -----------
                                  4           4
          (%o4)                 --------------------
                                         %pi 3/2
                                    (1 - ---)
                                          4

 -- Function: kurtosis_rayleigh (B)
     Returns the kurtosis coefficient of a Rayleigh(b) random variable, with b>0.

     The Rayleigh(b) random variable is equivalent to the Weibull(2,1/b), therefore when Maxima has not enough information to get the result, a noun form based on the Weibull kurtosis coefficient is returned.

          (%i1) load (distrib)$
          (%i2) kurtosis_rayleigh(b);
                                                   1
          (%o2)                kurtosis_weibull(2, -)
                                                   b
          (%i3) assume(b>0)$ kurtosis_rayleigh(b);
                                            2
                                       3 %pi
                                   2 - ------
                                         16
          (%o4)                    ---------- - 3
                                        %pi 2
                                   (1 - ---)
                                         4

 -- Function: random_rayleigh (B)
 -- Function: random_rayleigh (B,N)
     Returns a Rayleigh(b) random variate, with b>0. Calling random_rayleigh with a second argument N, a random sample of size N will be simulated.

     Only the inverse method is implemented. To make use of this function, write first load(distrib).

 -- Function: pdf_laplace (X,A,B)
     Returns the value at X of the density function of a Laplace(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: cdf_laplace (X,A,B)
     Returns the value at X of the distribution function of a Laplace(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: quantile_laplace (Q,A,B)
     Returns the Q-quantile of a Laplace(a,b) random variable, with b>0; in other words, this is the inverse of cdf_laplace. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_laplace (A,B)
     Returns the mean of a Laplace(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: var_laplace (A,B)
     Returns the variance of a Laplace(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: std_laplace (A,B)
     Returns the standard deviation of a Laplace(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: skewness_laplace (A,B)
     Returns the skewness coefficient of a Laplace(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: kurtosis_laplace (A,B)
     Returns the kurtosis coefficient of a Laplace(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: random_laplace (A,B)
 -- Function: random_laplace (A,B,N)
     Returns a Laplace(a,b) random variate, with b>0. Calling random_laplace with a third argument N, a random sample of size N will be simulated.

     Only the inverse method is implemented. To make use of this function, write first load(distrib).

 -- Function: pdf_cauchy (X,A,B)
     Returns the value at X of the density function of a Cauchy(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: cdf_cauchy (X,A,B)
     Returns the value at X of the distribution function of a Cauchy(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: quantile_cauchy (Q,A,B)
     Returns the Q-quantile of a Cauchy(a,b) random variable, with b>0; in other words, this is the inverse of cdf_cauchy. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: random_cauchy (A,B)
 -- Function: random_cauchy (A,B,N)
     Returns a Cauchy(a,b) random variate, with b>0. Calling random_cauchy with a third argument N, a random sample of size N will be simulated.

     Only the inverse method is implemented. To make use of this function, write first load(distrib).

 -- Function: pdf_gumbel (X,A,B)
     Returns the value at X of the density function of a Gumbel(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: cdf_gumbel (X,A,B)
     Returns the value at X of the distribution function of a Gumbel(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: quantile_gumbel (Q,A,B)
     Returns the Q-quantile of a Gumbel(a,b) random variable, with b>0; in other words, this is the inverse of cdf_gumbel. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_gumbel (A,B)
     Returns the mean of a Gumbel(a,b) random variable, with b>0.

          (%i1) load (distrib)$
          (%i2) assume(b>0)$  mean_gumbel(a,b);
          (%o3)                     %gamma b + a

     where symbol %gamma stands for the Euler-Mascheroni constant. See also %gamma.

 -- Function: var_gumbel (A,B)
     Returns the variance of a Gumbel(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: std_gumbel (A,B)
     Returns the standard deviation of a Gumbel(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: skewness_gumbel (A,B)
     Returns the skewness coefficient of a Gumbel(a,b) random variable, with b>0.

          (%i1) load (distrib)$
          (%i2) assume(b>0)$ skewness_gumbel(a,b);
                                 12 sqrt(6) zeta(3)
          (%o3)                  ------------------
                                           3
                                        %pi
          (%i4) numer:true$ skewness_gumbel(a,b);
          (%o5)                   1.139547099404649

     where zeta stands for the Riemann's zeta function.

 -- Function: kurtosis_gumbel (A,B)
     Returns the kurtosis coefficient of a Gumbel(a,b) random variable, with b>0. To make use of this function, write first load(distrib).

 -- Function: random_gumbel (A,B)
 -- Function: random_gumbel (A,B,N)
     Returns a Gumbel(a,b) random variate, with b>0. Calling random_gumbel with a third argument N, a random sample of size N will be simulated.

     Only the inverse method is implemented. To make use of this function, write first load(distrib).

Functions and Variables for discrete distributions
==================================================

 -- Function: pdf_binomial (X,N,P)
     Returns the value at X of the probability function of a Binomial(n,p) random variable, with 0<p<1 and n a positive integer. To make use of this function, write first load(distrib).

 -- Function: cdf_binomial (X,N,P)
     Returns the value at X of the distribution function of a Binomial(n,p) random variable, with 0<p<1 and n a positive integer.

     This function is numerically computed if the global variable numer equals true,  otherwise it returns a nominal expression.

          (%i1) load (distrib)$
          (%i2) cdf_binomial(5,7,1/6);
                                                   1
          (%o2)                 cdf_binomial(5, 7, -)
                                                   6
          (%i3) cdf_binomial(5,7,1/6), numer;
          (%o3)                   .9998713991769548

 -- Function: quantile_binomial (Q,N,P)
     Returns the Q-quantile of a Binomial(n,p) random variable, with 0<p<1 and n a positive integer; in other words, this is the inverse of cdf_binomial. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_binomial (N,P)
     Returns the mean of a Binomial(n,p) random variable, with 0<p<1 and n a positive integer. To make use of this function, write first load(distrib).

 -- Function: var_binomial (N,P)
     Returns the variance of a Binomial(n,p) random variable, with 0<p<1 and n a positive integer. To make use of this function, write first load(distrib).

 -- Function: std_binomial (N,P)
     Returns the standard deviation of a Binomial(n,p) random variable, with 0<p<1 and n a positive integer. To make use of this function, write first load(distrib).

 -- Function: skewness_binomial (N,P)
     Returns the skewness coefficient of a Binomial(n,p) random variable, with 0<p<1 and n a positive integer. To make use of this function, write first load(distrib).

 -- Function: kurtosis_binomial (N,P)
     Returns the kurtosis coefficient of a Binomial(n,p) random variable, with 0<p<1 and n a positive integer. To make use of this function, write first load(distrib).

 -- Option variable: random_binomial_algorithm
     Default value: kachit

     This is the selected algorithm for simulating random binomial variates. Implemented algorithms are kachit, bernoulli and inverse:

       * kachit, based on algorithm described in Kachitvichyanukul, V. and Schmeiser, B.W. (1988) BINOMIAL RANDOM VARIATE GENERATION. Communications of the ACM, 31, Feb., 216.

       * bernoulli, based on simulation of Bernoulli trials.

       * inverse, based on the general inverse method.

     See also random_binomial.

 -- Function: random_binomial (N,P)
 -- Function: random_binomial (N,P,M)
     Returns a Binomial(n,p) random variate, with 0<p<1 and n a positive integer. Calling random_binomial with a third argument M, a random sample of size M will be simulated.

     There are three algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_binomial_algorithm, which defaults to kachit.

     See also random_binomial_algorithm. To make use of this function, write first load(distrib).

 -- Function: pdf_poisson (X,M)
     Returns the value at X of the probability function of a Poisson(m) random variable, with m>0. To make use of this function, write first load(distrib).

 -- Function: cdf_poisson (X,M)
     Returns the value at X of the distribution function of a Poisson(m) random variable, with m>0.

     This function is numerically computed if the global variable numer equals true,  otherwise it returns a nominal expression.

          (%i1) load (distrib)$
          (%i2) cdf_poisson(3,5);
          (%o2)                   cdf_poisson(3, 5)
          (%i3) cdf_poisson(3,5), numer;
          (%o3)                   .2650259152973617

 -- Function: quantile_poisson (Q,M)
     Returns the Q-quantile of a Poisson(m) random variable, with m>0; in other words, this is the inverse of cdf_poisson. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_poisson (M)
     Returns the mean of a Poisson(m) random variable, with  m>0. To make use of this function, write first load(distrib).

 -- Function: var_poisson (M)
     Returns the variance of a Poisson(m) random variable, with  m>0. To make use of this function, write first load(distrib).

 -- Function: std_poisson (M)
     Returns the standard deviation of a Poisson(m) random variable, with m>0. To make use of this function, write first load(distrib).

 -- Function: skewness_poisson (M)
     Returns the skewness coefficient of a Poisson(m) random variable, with m>0. To make use of this function, write first load(distrib).

 -- Function: kurtosis_poisson (M)
     Returns the kurtosis coefficient of a Poisson random variable  Poi(m), with m>0. To make use of this function, write first load(distrib).

 -- Option variable: random_poisson_algorithm
     Default value: ahrens_dieter

     This is the selected algorithm for simulating random Poisson variates. Implemented algorithms are ahrens_dieter and inverse:

       * ahrens_dieter, based on algorithm described in Ahrens, J.H. and Dieter, U. (1982) COMPUTER GENERATION OF POISSON DEVIATES FROM MODIFIED NORMAL DISTRIBUTIONS. ACM Trans. Math. Software, 8, 2, June,163-179.

       * inverse, based on the general inverse method.

     See also random_poisson.

 -- Function: random_poisson (M)
 -- Function: random_poisson (M,N)
     Returns a Poisson(m) random variate, with m>0. Calling random_poisson with a second argument N, a random sample of size N will be simulated.

     There are two algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_poisson_algorithm, which defaults to ahrens_dieter.

     See also random_poisson_algorithm. To make use of this function, write first load(distrib).

 -- Function: pdf_bernoulli (X,P)
     Returns the value at X of the probability function of a Bernoulli(p) random variable, with 0<p<1.

     The Bernoulli(p) random variable is equivalent to the Binomial(1,p), therefore when Maxima has not enough information to get the result, a noun form based on the binomial probability function is returned.

          (%i1) load (distrib)$
          (%i2) pdf_bernoulli(1,p);
          (%o2)                 pdf_binomial(1, 1, p)
          (%i3) assume(0<p,p<1)$ pdf_bernoulli(1,p);
          (%o4)                           p

 -- Function: cdf_bernoulli (X,P)
     Returns the value at X of the distribution function of a Bernoulli(p) random variable, with 0<p<1. To make use of this function, write first load(distrib).

 -- Function: quantile_bernoulli (Q,P)
     Returns the Q-quantile of a Bernoulli(p) random variable, with 0<p<1; in other words, this is the inverse of cdf_bernoulli. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_bernoulli (P)
     Returns the mean of a Bernoulli(p) random variable, with 0<p<1.

     The Bernoulli(p) random variable is equivalent to the Binomial(1,p), therefore when Maxima has not enough information to get the result, a noun form based on the binomial mean is returned.

          (%i1) load (distrib)$
          (%i2) mean_bernoulli(p);
          (%o2)                  mean_binomial(1, p)
          (%i3) assume(0<p,p<1)$ mean_bernoulli(p);
          (%o4)                           p

 -- Function: var_bernoulli (P)
     Returns the variance of a Bernoulli(p) random variable, with 0<p<1.

     The Bernoulli(p) random variable is equivalent to the Binomial(1,p), therefore when Maxima has not enough information to get the result, a noun form based on the binomial variance is returned.

          (%i1) load (distrib)$
          (%i2) var_bernoulli(p);
          (%o2)                  var_binomial(1, p)
          (%i3) assume(0<p,p<1)$ var_bernoulli(p);
          (%o4)                       (1 - p) p

 -- Function: std_bernoulli (P)
     Returns the standard deviation of a Bernoulli(p) random variable, with 0<p<1.

     The Bernoulli(p) random variable is equivalent to the Binomial(1,p), therefore when Maxima has not enough information to get the result, a noun form based on the binomial standard deviation is returned.

          (%i1) load (distrib)$
          (%i2) std_bernoulli(p);
          (%o2)                  std_binomial(1, p)
          (%i3) assume(0<p,p<1)$ std_bernoulli(p);
          (%o4)                  sqrt(1 - p) sqrt(p)

 -- Function: skewness_bernoulli (P)
     Returns the skewness coefficient of a Bernoulli(p) random variable, with 0<p<1.

     The Bernoulli(p) random variable is equivalent to the Binomial(1,p), therefore when Maxima has not enough information to get the result, a noun form based on the binomial skewness coefficient is returned.

          (%i1) load (distrib)$
          (%i2) skewness_bernoulli(p);
          (%o2)                skewness_binomial(1, p)
          (%i3) assume(0<p,p<1)$ skewness_bernoulli(p);
                                       1 - 2 p
          (%o4)                  -------------------
                                 sqrt(1 - p) sqrt(p)

 -- Function: kurtosis_bernoulli (P)
     Returns the kurtosis coefficient of a Bernoulli(p) random variable, with 0<p<1.

     The Bernoulli(p) random variable is equivalent to the Binomial(1,p), therefore when Maxima has not enough information to get the result, a noun form based on the binomial kurtosis coefficient is returned.

          (%i1) load (distrib)$
          (%i2) kurtosis_bernoulli(p);
          (%o2)                kurtosis_binomial(1, p)
          (%i3) assume(0<p,p<1)$ kurtosis_bernoulli(p);
                                   1 - 6 (1 - p) p
          (%o4)                    ---------------
                                      (1 - p) p

 -- Function: random_bernoulli (P)
 -- Function: random_bernoulli (P,N)
     Returns a Bernoulli(p) random variate, with 0<p<1. Calling random_bernoulli with a second argument N, a random sample of size N will be simulated.

     This is a direct application of the random built-in Maxima function.

     See also random. To make use of this function, write first load(distrib).

 -- Function: pdf_geometric (X,P)
     Returns the value at X of the probability function of a Geometric(p) random variable, with 0<p<1. To make use of this function, write first load(distrib).

 -- Function: cdf_geometric (X,P)
     Returns the value at X of the distribution function of a Geometric(p) random variable, with 0<p<1. To make use of this function, write first load(distrib).

 -- Function: quantile_geometric (Q,P)
     Returns the Q-quantile of a Geometric(p) random variable, with 0<p<1; in other words, this is the inverse of cdf_geometric. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_geometric (P)
     Returns the mean of a Geometric(p) random variable, with 0<p<1. To make use of this function, write first load(distrib).

 -- Function: var_geometric (P)
     Returns the variance of a Geometric(p) random variable, with 0<p<1. To make use of this function, write first load(distrib).

 -- Function: std_geometric (P)
     Returns the standard deviation of a Geometric(p) random variable, with 0<p<1. To make use of this function, write first load(distrib).

 -- Function: skewness_geometric (P)
     Returns the skewness coefficient of a Geometric(p) random variable, with 0<p<1. To make use of this function, write first load(distrib).

 -- Function: kurtosis_geometric (P)
     Returns the kurtosis coefficient of a geometric random variable  Geo(p), with 0<p<1. To make use of this function, write first load(distrib).

 -- Option variable: random_geometric_algorithm
     Default value: bernoulli

     This is the selected algorithm for simulating random geometric variates. Implemented algorithms are bernoulli, devroye and inverse:

       * bernoulli, based on simulation of Bernoulli trials.

       * devroye, based on algorithm described in Devroye, L. (1986) NON-UNIFORM RANDOM VARIATE GENERATION. Springer Verlag, p. 480.

       * inverse, based on the general inverse method.

     See also random_geometric.

 -- Function: random_geometric (P)
 -- Function: random_geometric (P,N)
     Returns a Geometric(p) random variate, with 0<p<1. Calling random_geometric with a second argument N, a random sample of size N will be simulated.

     There are three algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_geometric_algorithm, which defaults to bernoulli.

     See also random_geometric_algorithm. To make use of this function, write first load(distrib).

 -- Function: pdf_discrete_uniform (X,N)
     Returns the value at X of the probability function of a Discrete Uniform(n) random variable, with n a strictly positive integer. To make use of this function, write first load(distrib).

 -- Function: cdf_discrete_uniform (X,N)
     Returns the value at X of the distribution function of a Discrete Uniform(n) random variable, with n a strictly positive integer. To make use of this function, write first load(distrib).

 -- Function: quantile_discrete_uniform (Q,N)
     Returns the Q-quantile of a Discrete Uniform(n) random variable, with n a strictly positive integer; in other words, this is the inverse of cdf_discrete_uniform. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_discrete_uniform (N)
     Returns the mean of a Discrete Uniform(n) random variable, with n a strictly positive integer. To make use of this function, write first load(distrib).

 -- Function: var_discrete_uniform (N)
     Returns the variance of a Discrete Uniform(n) random variable, with n a strictly positive integer. To make use of this function, write first load(distrib).

 -- Function: std_discrete_uniform (N)
     Returns the standard deviation of a Discrete Uniform(n) random variable, with n a strictly positive integer. To make use of this function, write first load(distrib).

 -- Function: skewness_discrete_uniform (N)
     Returns the skewness coefficient of a Discrete Uniform(n) random variable, with n a strictly positive integer. To make use of this function, write first load(distrib).

 -- Function: kurtosis_discrete_uniform (N)
     Returns the kurtosis coefficient of a Discrete Uniform(n) random variable, with n a strictly positive integer. To make use of this function, write first load(distrib).

 -- Function: random_discrete_uniform (N)
 -- Function: random_discrete_uniform (N,M)
     Returns a Discrete Uniform(n) random variate, with n a strictly positive integer. Calling random_discrete_uniform with a second argument M, a random sample of size M will be simulated.

     This is a direct application of the random built-in Maxima function.

     See also random. To make use of this function, write first load(distrib).

 -- Function: pdf_hypergeometric (X,N1,N2,N)
     Returns the value at X of the probability function of a Hypergeometric(n1,n2,n) random variable, with N1, N2 and N non negative integers and n<=n1+n2. To make use of this function, write first load(distrib).

 -- Function: cdf_hypergeometric (X,N1,N2,N)
     Returns the value at X of the distribution function of a Hypergeometric(n1,n2,n) random variable, with N1, N2 and N non negative integers and n<=n1+n2. To make use of this function, write first load(distrib).

 -- Function: quantile_hypergeometric (Q,N1,N2,N)
     Returns the Q-quantile of a Hypergeometric(n1,n2,n) random variable, with N1, N2 and N non negative integers and n<=n1+n2; in other words, this is the inverse of cdf_hypergeometric. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_hypergeometric (N1,N2,N)
     Returns the mean of a discrete uniform random variable Hyp(n1,n2,n), with N1, N2 and N non negative integers and n<=n1+n2. To make use of this function, write first load(distrib).

 -- Function: var_hypergeometric (N1,N2,N)
     Returns the variance of a hypergeometric  random variable Hyp(n1,n2,n), with N1, N2 and N non negative integers and n<=n1+n2. To make use of this function, write first load(distrib).

 -- Function: std_hypergeometric (N1,N2,N)
     Returns the standard deviation of a Hypergeometric(n1,n2,n) random variable, with N1, N2 and N non negative integers and n<=n1+n2. To make use of this function, write first load(distrib).

 -- Function: skewness_hypergeometric (N1,N2,N)
     Returns the skewness coefficient of a Hypergeometric(n1,n2,n) random variable, with N1, N2 and N non negative integers and n<=n1+n2. To make use of this function, write first load(distrib).

 -- Function: kurtosis_hypergeometric (N1,N2,N)
     Returns the kurtosis coefficient of a Hypergeometric(n1,n2,n) random variable, with N1, N2 and N non negative integers and n<=n1+n2. To make use of this function, write first load(distrib).

 -- Option variable: random_hypergeometric_algorithm
     Default value: kachit

     This is the selected algorithm for simulating random hypergeometric variates. Implemented algorithms are kachit and inverse:

       * kachit, based on algorithm described in Kachitvichyanukul, V., Schmeiser, B.W. (1985) COMPUTER GENERATION OF HYPERGEOMETRIC RANDOM VARIATES. Journal of Statistical Computation and Simulation 22, 127-145.

       * inverse, based on the general inverse method.

     See also random_hypergeometric.

 -- Function: random_hypergeometric (N1,N2,N)
 -- Function: random_hypergeometric (N1,N2,N,M)
     Returns a Hypergeometric(n1,n2,n) random variate, with N1, N2 and N non negative integers and n<=n1+n2. Calling random_hypergeometric with a fourth argument M, a random sample of size M will be simulated.

     There are two algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_hypergeometric_algorithm, which defaults to kachit.

     See also random_hypergeometric_algorithm. To make use of this function, write first load(distrib).

 -- Function: pdf_negative_binomial (X,N,P)
     Returns the value at X of the probability function of a Negative Binomial(n,p) random variable, with 0<p<1 and n a positive integer. To make use of this function, write first load(distrib).

 -- Function: cdf_negative_binomial (X,N,P)
     Returns the value at X of the distribution function of a Negative Binomial(n,p) random variable, with 0<p<1 and n a positive integer.

     This function is numerically computed if the global variable numer equals true,  otherwise it returns a nominal expression.

          (%i1) load (distrib)$
          (%i2) cdf_negative_binomial(3,4,1/8);
                                                       1
          (%o2)            cdf_negative_binomial(3, 4, -)
                                                       8
          (%i3) cdf_negative_binomial(3,4,1/8), numer;
          (%o3)                  .006238937377929698

 -- Function: quantile_negative_binomial (Q,N,P)
     Returns the Q-quantile of a Negative Binomial(n,p) random variable, with 0<p<1 and n a positive integer; in other words, this is the inverse of cdf_negative_binomial. Argument Q must be an element of [0,1]. To make use of this function, write first load(distrib).

 -- Function: mean_negative_binomial (N,P)
     Returns the mean of a Negative Binomial(n,p) random variable, with 0<p<1 and n a positive integer. To make use of this function, write first load(distrib).

 -- Function: var_negative_binomial (N,P)
     Returns the variance of a Negative Binomial(n,p) random variable, with 0<p<1 and n a positive integer. To make use of this function, write first load(distrib).

 -- Function: std_negative_binomial (N,P)
     Returns the standard deviation of a Negative Binomial(n,p) random variable, with 0<p<1 and n a positive integer. To make use of this function, write first load(distrib).

 -- Function: skewness_negative_binomial (N,P)
     Returns the skewness coefficient of a Negative Binomial(n,p) random variable, with 0<p<1 and n a positive integer. To make use of this function, write first load(distrib).

 -- Function: kurtosis_negative_binomial (N,P)
     Returns the kurtosis coefficient of a Negative Binomial(n,p) random variable, with 0<p<1 and n a positive integer. To make use of this function, write first load(distrib).

 -- Option variable: random_negative_binomial_algorithm
     Default value: bernoulli

     This is the selected algorithm for simulating random negative binomial variates. Implemented algorithms are devroye, bernoulli and inverse:

       * devroye, based on algorithm described in Devroye, L. (1986) NON-UNIFORM RANDOM VARIATE GENERATION. Springer Verlag, p. 480.

       * bernoulli, based on simulation of Bernoulli trials.

       * inverse, based on the general inverse method.

     See also random_negative_binomial.

 -- Function: random_negative_binomial (N,P)
 -- Function: random_negative_binomial (N,P,M)
     Returns a Negative Binomial(n,p) random variate, with 0<p<1 and n a positive integer. Calling random_negative_binomial with a third argument M, a random sample of size M will be simulated.

     There are three algorithms implemented for this function, the one to be used can be selected giving a certain value to the global variable random_negative_binomial_algorithm, which defaults to bernoulli.

     See also random_negative_binomial_algorithm. To make use of this function, write first load(distrib).