F – F-distribution
beta – Beta distribution
cauchy – Cauchy distribution
chi – Chi distribution
chisquare – Chisquare distribution
exponential – Exponential distribution
extremeI – Extreme value type I (Gumbel-type) distribution
extremeII – Extreme value type II (Frechet-type) distribution
gamma – Gamma distribution
gig – Generalized Inverse Gaussian distribution
gig2 – Generalized Inverse Gaussian distribution
hyperbolic – Hyperbolic distribution
ig – Inverse Gaussian distribution
laplace – Laplace distribution
logistic – Logistic distribution
lognormal – Log-Normal distribution
lomax – Lomax distribution (Pareto distribution of second kind)
normal – Normal distribution
pareto – Pareto distribution (of first kind)
powerexponential – Powerexponential (Subbotin) distribution
rayleigh – Rayleigh distribution
slash – Slash distribution
student – Student's t distribution
triangular – Triangular distribution
uniform – Uniform distribution
weibull – Weibull distribution
UNU.RAN (Universal Non-Uniform RAndom Number generator) is a collection of algorithms for generating non-uniform pseudorandom variates as a library of C functions designed and implemented by the ARVAG (Automatic Random VAriate Generation) project group in Vienna, and released under the GNU Public License (GPL). It is especially designed for such situations where
Of course it is also well suited for standard distributions. However due to its more sophisticated programming interface it might not be as easy to use if you only look for a generator for the standard normal distribution. (Although UNU.RAN provides generators that are superior in many aspects to those found in quite a number of other libraries.)
UNU.RAN implements several methods for generating random numbers. The choice depends primary on the information about the distribution can be provided and – if the user is familar with the different methods – on the preferences of the user.
The design goals of UNU.RAN are to provide reliable, portable and robust (as far as this is possible) functions with a consisent and easy to use interface. It is suitable for all situation where experiments with different distributions including non-standard distributions. For example it is no problem to replace the normal distribution by an empirical distribution in a model.
Since originally designed as a library for so called black-box or universal algorithms its interface is different from other libraries. (Nevertheless it also contains special generators for standard distributions.) It does not provide subroutines for random variate generation for particular distributions. Instead it uses an object-oriented interface. Distributions and generators are treated as independent objects. This approach allows one not only to have different methods for generating non-uniform random variates. It is also possible to choose the method which is optimal for a given situation (e.g. speed, quality of random numbers, using for variance reduction techniques, etc.). It also allows to sample from non-standard distribution or even from distributions that arise in a model and can only be computed in a complicated subroutine.
Sampling from a particular distribution requires the following steps:
There are four types of objects that can be manipulated independently:
Of course a library of standard distributions is included (and these can be further modified to get, e.g., truncated distributions). Moreover the library provides subroutines to build almost arbitrary distributions.
NULL pointer is returned in the initialization step).
We designed this document in a way such that one can use UNU.RAN with reading as little as necessary. Read Installation for the instructions to install the library. Concepts of UNU.RAN, discribes the basics of UNU.RAN. It also has a short guideline for choosing an appropriate method. In Examples examples are given that can be copied and modified. They also can be found in the directory examples in the source tree.
Further information are given in consecutive chapters. Handling distribution objects, describes how to create and manipulate distribution objects. standard distributions, describes predefined distribution objects that are ready to use. Methods describes the various methods in detail. For each of possible distribution classes (continuous, discrete, empirical, multivariate) there exists a short overview section that can be used to choose an appropriate method followed by sections that describe each of the particular methods in detail. These are merely for users with some knowledge about the methods who want to change method-specific parameters and can be ignored by others.
Abbreviations and explanation of some basic terms can be found in Glossary.
UNU.RAN was developed on an Intel architecture under Linux with the GNU C compiler but should compile and run on any computing environment. It requires an ANSI compliant C compiler.
Below find the installation instructions for unices.
UNU.RAN can be used with any uniform random number generator but (at the moment) some features work best with Pierre L'Ecuyer's RngStreams library (see http://statmath.wu.ac.at/software/RngStreams/ for a description and downloading. For details on using uniform random number in UNU.RAN see Using uniform random number generators.
Install the required libraries first.
tar zxvf unuran-1.6.0.tar.gz
cd unuran-1.6.0
sh ./configure --prefix=<prefix>
where <prefix> is the root of the installation tree.
When omitted /usr/local is used.
Use ./configure --help to get a list of other options.
In particular the following flags are important:
--with-urng-rngstreamno]
--with-urng-prngno]
--with-urng-gslno]
--with-urng-defaultbuiltin]
We strongly recommend to use RngStreams library:
sh ./configure --with-urng-rngstream --with-urng-default=rngstream
Important: You must install the respective libraries
RngStreams, PRNG and GSL before
./configure is executed.
--enable-sharedno]
unur_gen_info for
information about generator objects. This is intented for
using in interactive computing environments.
This feature can be enabled / disabled by means of the
configure flag
--enable-infoyes]
--enable-deprecatedno]
--enable-check-structno]
--enable-logging make
make install
Obviously $(prefix)/include and $(prefix)/lib
must be in the search path of your compiler. You can use environment
variables to add these directories to the search path. If you
are using the bash type (or add to your profile):
export LIBRARY_PATH="<prefix>/lib"
export C_INCLURE_PATH="<prefix>/include"
If you want to make a shared library, then making such a library can be enabled using
sh ./configure --enable-shared
If you want to link against the shared library make sure that
it can be found when executing the binary that links to the
library. If it is not installed in the usual path, then the
easiest way is to set the LD_LIBRARY_PATH environment
variable. See any operating system documentation about shared
libraries for more information, such as the ld(1) and
ld.so(8) manual pages.
make check
However, some of these tests requires the usage of the PRNG or RngStreams library and are only executed if these are installed enabled by the corresponding configure flag.
An extended set of tests is run by
make fullcheck
However some of these might fail occasionally due to roundoff errors or the mysteries of floating point arithmetic, since we have used some extreme settings to test the library.
UNU.RAN now relies on some aspects of IEEE
754 compliant floating point arithmetic. In particular,
1./0. and 0./0. must result in infinity
and NaN (not a number), respectively, and must not
cause a floating point exception.
For allmost all modern compting architecture this is implemented
in hardware. For others there should be a special compiler flag
to get this feature (e.g., -MIEEE on DEC alpha or
-mp for the Intel C complier).
With UNU.RAN version 1.0.x some of the macro definitions in
file src/unuran_config.h are moved into file
config.h and are set/controlled by the
./configure script.
Writting logging information into the logfile must now be enabled when running the configure script:
sh ./configure --enable-logging
With UNU.RAN version 0.8.0 the interface for changing underlying distributions and running a reinitialization routine has been simplified. The old routines can be compiled into the library using the following configure flag:
sh ./configure --enable-deprecated
Notice: Using these deprecated routines is not supported any more and this strong discouraged.
Wrapper functions for external sources of uniform random numbers are now enabled by configure flags and not by macros defined in file src/unuran_config.h.
The file src/unuran_config.h is not installed any more. It is now only included when the library is compiled. It should be removed from the global include path of the compiler.
The library is written in ANSI C and is intended to conform to the ANSI C standard. It should be portable to any system with a working ANSI C compiler.
The library does not rely on any non-ANSI extensions in the interface it exports to the user. Programs you write using UNU.RAN can be ANSI compliant. Extensions which can be used in a way compatible with pure ANSI C are supported, however, via conditional compilation. This allows the library to take advantage of compiler extensions on those platforms which support them.
To avoid namespace conflicts all exported function names and
variables have the prefix unur_, while exported macros have
the prefix UNUR_.
If you want to use the library you must include the UNU.RAN header file
#include <unuran.h>
If you also need the test routines then also add
#include <unuran_tests.h>
If wrapper functions for external sources of uniform random number generators are used, the corresponding header files must also be included, e.g.,
#include <unuran_urng_rngstream.h>
If these header files are not installed on the standard search path
of your compiler you will also need to provide its location to the
preprocessor as a command line flag. The default location of the
unuran.h is /usr/local/include. A typical compilation
command for a source file app.c with the GNU C compiler
gcc is,
gcc -I/usr/local/include -c app.c
This results in an object file app.o. The default include
path for gcc searches /usr/local/include
automatically so the -I option can be omitted when UNU.RAN is
installed in its default location.
The library is installed as a single file, libunuran.a. A shared version of the library is also installed on systems that support shared libraries. The default location of these files is /usr/local/lib. To link against the library you need to specify the main library. The following example shows how to link an application with the library (and the the RNGSTREAMS library if you decide to use this source of uniform pseudo-random numbers),
gcc app.o -lunuran -lrngstreams -lm
To run a program linked with the shared version of the library it
may be necessary to define the shell variable
LD_LIBRARY_PATH to include the directory where the library
is installed. For example,
LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH
To compile a statically linked version of the program instead, use the
-static flag in gcc,
gcc -static app.o -lunuran -lrngstreams -lm
The library header files automatically define functions to have
extern "C" linkage when included in C++ programs.
UNU.RAN is a C library for generating non-uniformly distributed random variates. Its emphasis is on the generation of non-standard distribution and on streams of random variates of special purposes. It is designed to provide a consistent tool to sample from distributions with various properties. Since there is no universal method that fits for all situations, various methods for sampling are implemented.
UNU.RAN solves this complex task by means of an object oriented programming interface. Three basic objects are used:
UNUR_DISTRUNUR_GENUNUR_PARNotice that the parameter objects only hold pointers to arrays but do not have their own copy of such an array. Especially, if a dynamically allocated array is used it must not be freed until the generator object has been created!
The idea behind these structures is that creatin distributions, choosing a generation method and draing samples are orthogonal (ie. independent) functions of the library. The parameter object is only introduced due to the necessity to deal with various parameters and switches for each of these generation methods which are required to adjust the algorithms to unusual distributions with extreme properties but have default values that are suitable for most applications. These parameters and the data for distributions are set by various functions.
Once a generator object has been created sampling (from the univariate continuous distribution) can be done by the following command:
double x = unur_sample_cont(generator);
Analogous commands exist for discrete and multivariate distributions. For detailed examples that can be copied and modified see Examples.
All information about a distribution are stored in objects
(structures) of type UNUR_DISTR.
UNU.RAN has five different types of distribution objects:
contcvecdiscrcempcvempmatrDistribution objects can be created from scratch by the following call
distr = unur_distr_<type>_new();
where <type> is one of the five possible types from the
above table.
Notice that these commands only create an empty object which
still must be filled by means of calls for each type of
distribution object
(see Handling distribution objects).
The naming scheme of these functions is designed to indicate the
corresponding type of the distribution object and the task to be
performed. It is demonstated on the following example.
unur_distr_cont_set_pdf(distr, mypdf);
This command stores a PDF named mypdf in the distribution
object distr which must have the type cont.
Of course UNU.RAN provides an easier way to use standard distributions. Instead of using unur_distr_<type>_new calls and fuctions unur_distr_<type>_set_<...> for setting data, objects for standard distribution can be created by a single call. Eg. to get an object for the normal distribution with mean 2 and standard deviation 5 use
double parameter[2] = {2.0 ,5.0};
UNUR_DISTR *distr = unur_distr_normal(parameter, 2);
For a list of standard distributions see Standard distributions.
The information that a distribution object must contain depends heavily on the chosen generation method choosen.
Brackets indicate optional information while a tilde indicates that only an approximation must be provided. See Glossary, for unfamiliar terms.
|
Methods for continuous univariate distributions sample with unur_sample_cont
|
|
Methods for continuous empirical univariate distributions sample with unur_sample_cont EMPK: Requires an observed sample. EMPL: Requires an observed sample. |
|
Methods for continuous multivariate distributions sample with unur_sample_vec NORTA: Requires rank correlation matrix and marginal distributions. VNROU: Requires the PDF. MVSTD: Generator for built-in standard distributions. MVTDR: Requires PDF and gradiant of PDF. |
|
Methods for continuous empirical multivariate distributions sample with unur_sample_vec VEMPK: Requires an observed sample. |
|
Methods for discrete univariate distributions sample with unur_sample_discr
|
|
Methods for matrix distributions sample with unur_sample_matr MCORR: Distribution object for random correlation matrix. |
|
Markov Chain Methods for continuous multivariate distributions sample with unur_sample_vec GIBBS: T-concave logPDF and derivatives of logPDF. HITRO: Requires PDF. |
Because of tremendous variety of possible problems, UNU.RAN provides many methods. All information for creating a generator object has to be collected in a parameter object first. For example, if the task is to sample from a continuous distribution the method AROU might be a good choice. Then the call
UNUR_PAR *par = unur_arou_new(distribution);
creates an parameter object par with a pointer to the
distribution object and default values for all necessary parameters
for method AROU.
Other methods can be used by replacing arou with the name
of the desired methods (in lower case letters):
UNUR_PAR *par = unur_<method>_new(distribution);
This sets the default values for all necessary parameters for the chosen method. These are suitable for almost all applications. Nevertheless, it is possible to control the behavior of the method using corresponding set calls for each method. This might be necessary to adjust the algorithm for an unusual distribution with extreme properties, or just for fine tuning the perforence of the algorithm. The following example demonstrates how to change the maximum number of iterations for method NINV to the value 50:
unur_ninv_set_max_iteration(par, 50);
All available methods are described in details in Methods.
Now it is possible to create a generator object:
UNUR_GEN *generator = unur_init(par);
if (generator == NULL) exit(EXIT_FAILURE);
Important: You must always check whether
unur_init
has
been executed successfully. Otherwise the NULL pointer is returned
which causes a segmentation fault when used for sampling.
Important:
The call of
unur_init
destroys the parameter object!
Moreover, it is recommended to call
unur_init
immediately after
the parameter object par has created and modified.
An existing generator object is a rather static construct. Nevertheless, some of the parameters can still be modified by chg calls, e.g.
unur_ninv_chg_max_iteration(gen, 30);
Notice that it is important when parameters are changed because different functions must be used:
The function name includes the term set and the first
argument must be of type UNUR_PAR when the parameters are
changed before the generator object is created.
The function name includes the term chg and the first
argument must be of type UNUR_GEN when the parameters are
changed for an existing generator object.
For details see Methods.
You can now use your generator object in any place of your program
to sample from your distribution. You only have to take care about
the type of variates it computes: double, int or a
vector (array of doubles).
Notice that at this point it does not matter whether you are
sampling from a gamma distribution, a truncated normal distribution
or even an empirical distribution.
It is possible for a generator object to change the parameters and the domain of the underlying distribution. This must be done by extracting this object by means of a unur_get_distr call and changing the distribution using the correspondig set calls, see Handling distribution objects. The generator object must then be reinitialized by means of the unur_reinit call.
Important: Currently not all methods allow reinitialization, see the description of the particular method (keyword Reinit).
When you do not need your generator object any more, you should destroy it:
unur_free(generator);
Each generator object can have its own uniform random number generator or share one with others. When created a parameter object the pointer for the uniform random number generator is set to the default generator. However, it can be changed at any time to any other generator:
unur_set_urng(par, urng);
or
unur_chg_urng(generator, urng);
respectively. See Using uniform random number generators, for details.
If you have any problems with UNU.RAN, suggestions how to improve the library, or find a bug, please contact us via email unuran@statmath.wu.ac.at.
For news please visit out homepage at http://statmath.wu.ac.at/unuran/.
The examples in this chapter should compile cleanly and can be found in the directory examples of the source tree of UNU.RAN. Assuming that UNU.RAN as well as the PRNG libraries have been installed properly (see Installation) each of these can be compiled (using the GCC in this example) with
gcc -Wall -O2 -o example example.c -lunuran -lprng -lm
Remark: -lprng must be omitted when the PRNG library
is not installed. Then however some of the examples might not work.
The library uses three objects:
UNUR_DISTR, UNUR_PAR and UNUR_GEN.
It is not important to understand the details of these objects but
it is important not to changed the order of their creation.
The distribution object can be destroyed after the generator
object has been made. (The parameter object is freed automatically
by the
unur_init
call.) It is also important to check the result
of the
unur_init
call. If it has failed the NULL pointer is
returned and causes a segmentation fault when used for sampling.
We give all examples with the UNU.RAN standard API and the more convenient string API.
Select a distribution and let UNU.RAN do all necessary steps.
/* ------------------------------------------------------------- */
/* File: example0.c */
/* ------------------------------------------------------------- */
/* Include UNURAN header file. */
#include <unuran.h>
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double x; /* will hold the random number */
/* Declare the three UNURAN objects. */
UNUR_DISTR *distr; /* distribution object */
UNUR_PAR *par; /* parameter object */
UNUR_GEN *gen; /* generator object */
/* Use a predefined standard distribution: */
/* Gaussian with mean zero and standard deviation 1. */
/* Since this is the standard form of the distribution, */
/* we need not give these parameters. */
distr = unur_distr_normal(NULL, 0);
/* Use method AUTO: */
/* Let UNURAN select a suitable method for you. */
par = unur_auto_new(distr);
/* Now you can change some of the default settings for the */
/* parameters of the chosen method. We don't do it here. */
/* Create the generator object. */
gen = unur_init(par);
/* Notice that this call has also destroyed the parameter */
/* object `par' as a side effect. */
/* It is important to check if the creation of the generator */
/* object was successful. Otherwise `gen' is the NULL pointer */
/* and would cause a segmentation fault if used for sampling. */
if (gen == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* It is possible to reuse the distribution object to create */
/* another generator object. If you do not need it any more, */
/* it should be destroyed to free memory. */
unur_distr_free(distr);
/* Now you can use the generator object `gen' to sample from */
/* the standard Gaussian distribution. */
/* Eg.: */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* When you do not need the generator object any more, you */
/* can destroy it. */
unur_free(gen);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
Select a distribution and let UNU.RAN do all necessary steps.
/* ------------------------------------------------------------- */
/* File: example0_str.c */
/* ------------------------------------------------------------- */
/* String API. */
/* ------------------------------------------------------------- */
/* Include UNURAN header file. */
#include <unuran.h>
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double x; /* will hold the random number */
/* Declare UNURAN generator object. */
UNUR_GEN *gen; /* generator object */
/* Create the generator object. */
/* Use a predefined standard distribution: */
/* Standard Gaussian distribution. */
/* Use method AUTO: */
/* Let UNURAN select a suitable method for you. */
gen = unur_str2gen("normal()");
/* It is important to check if the creation of the generator */
/* object was successful. Otherwise `gen' is the NULL pointer */
/* and would cause a segmentation fault if used for sampling. */
if (gen == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* Now you can use the generator object `gen' to sample from */
/* the standard Gaussian distribution. */
/* Eg.: */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* When you do not need the generator object any more, you */
/* can destroy it. */
unur_free(gen);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
Select method AROU and use it with default parameters.
/* ------------------------------------------------------------- */
/* File: example1.c */
/* ------------------------------------------------------------- */
/* Include UNURAN header file. */
#include <unuran.h>
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double x; /* will hold the random number */
/* Declare the three UNURAN objects. */
UNUR_DISTR *distr; /* distribution object */
UNUR_PAR *par; /* parameter object */
UNUR_GEN *gen; /* generator object */
/* Use a predefined standard distribution: */
/* Gaussian with mean zero and standard deviation 1. */
/* Since this is the standard form of the distribution, */
/* we need not give these parameters. */
distr = unur_distr_normal(NULL, 0);
/* Choose a method: AROU. */
/* For other (suitable) methods replace "arou" with the */
/* respective name (in lower case letters). */
par = unur_arou_new(distr);
/* Now you can change some of the default settings for the */
/* parameters of the chosen method. We don't do it here. */
/* Create the generator object. */
gen = unur_init(par);
/* Notice that this call has also destroyed the parameter */
/* object `par' as a side effect. */
/* It is important to check if the creation of the generator */
/* object was successful. Otherwise `gen' is the NULL pointer */
/* and would cause a segmentation fault if used for sampling. */
if (gen == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* It is possible to reuse the distribution object to create */
/* another generator object. If you do not need it any more, */
/* it should be destroyed to free memory. */
unur_distr_free(distr);
/* Now you can use the generator object `gen' to sample from */
/* the standard Gaussian distribution. */
/* Eg.: */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* When you do not need the generator object any more, you */
/* can destroy it. */
unur_free(gen);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
Select method AROU and use it with default parameters.
/* ------------------------------------------------------------- */
/* File: example1_str.c */
/* ------------------------------------------------------------- */
/* String API. */
/* ------------------------------------------------------------- */
/* Include UNURAN header file. */
#include <unuran.h>
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double x; /* will hold the random number */
/* Declare UNURAN generator object. */
UNUR_GEN *gen; /* generator object */
/* Create the generator object. */
/* Use a predefined standard distribution: */
/* Standard Gaussian distribution. */
/* Choose a method: AROU. */
/* For other (suitable) methods replace "arou" with the */
/* respective name. */
gen = unur_str2gen("normal() & method=arou");
/* It is important to check if the creation of the generator */
/* object was successful. Otherwise `gen' is the NULL pointer */
/* and would cause a segmentation fault if used for sampling. */
if (gen == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* Now you can use the generator object `gen' to sample from */
/* the standard Gaussian distribution. */
/* Eg.: */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* When you do not need the generator object any more, you */
/* can destroy it. */
unur_free(gen);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
If you want to sample from a non-standard distribution, UNU.RAN might be exactly what you need. Depending on the information is available, a method must be choosen for sampling, see Concepts for an overview and Methods for details.
/* ------------------------------------------------------------- */
/* File: example2.c */
/* ------------------------------------------------------------- */
/* Include UNURAN header file. */
#include <unuran.h>
/* ------------------------------------------------------------- */
/* In this example we build a distribution object from scratch */
/* and sample from this distribution. */
/* */
/* We use method TDR (Transformed Density Rejection) which */
/* required a PDF and the derivative of the PDF. */
/* ------------------------------------------------------------- */
/* Define the PDF and dPDF of our distribution. */
/* */
/* Our distribution has the PDF */
/* */
/* / 1 - x*x if |x| <= 1 */
/* f(x) = < */
/* \ 0 otherwise */
/* */
/* The PDF of our distribution: */
double mypdf( double x, const UNUR_DISTR *distr )
/* The second argument (`distr') can be used for parameters */
/* for the PDF. (We do not use parameters in our example.) */
{
if (fabs(x) >= 1.)
return 0.;
else
return (1.-x*x);
} /* end of mypdf() */
/* The derivative of the PDF of our distribution: */
double mydpdf( double x, const UNUR_DISTR *distr )
{
if (fabs(x) >= 1.)
return 0.;
else
return (-2.*x);
} /* end of mydpdf() */
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double x; /* will hold the random number */
/* Declare the three UNURAN objects. */
UNUR_DISTR *distr; /* distribution object */
UNUR_PAR *par; /* parameter object */
UNUR_GEN *gen; /* generator object */
/* Create a new distribution object from scratch. */
/* It is a continuous distribution, and we need a PDF and the */
/* derivative of the PDF. Moreover we set the domain. */
/* Get empty distribution object for a continuous distribution */
distr = unur_distr_cont_new();
/* Assign the PDF and dPDF (defined above). */
unur_distr_cont_set_pdf( distr, mypdf );
unur_distr_cont_set_dpdf( distr, mydpdf );
/* Set the domain of the distribution (optional for TDR). */
unur_distr_cont_set_domain( distr, -1., 1. );
/* Choose a method: TDR. */
par = unur_tdr_new(distr);
/* Now you can change some of the default settings for the */
/* parameters of the chosen method. We don't do it here. */
/* Create the generator object. */
gen = unur_init(par);
/* Notice that this call has also destroyed the parameter */
/* object `par' as a side effect. */
/* It is important to check if the creation of the generator */
/* object was successful. Otherwise `gen' is the NULL pointer */
/* and would cause a segmentation fault if used for sampling. */
if (gen == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* It is possible to reuse the distribution object to create */
/* another generator object. If you do not need it any more, */
/* it should be destroyed to free memory. */
unur_distr_free(distr);
/* Now you can use the generator object `gen' to sample from */
/* the distribution. Eg.: */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* When you do not need the generator object any more, you */
/* can destroy it. */
unur_free(gen);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
If you want to sample from a non-standard distribution, UNU.RAN might be exactly what you need. Depending on the information is available, a method must be choosen for sampling, see Concepts for an overview and Methods for details.
/* ------------------------------------------------------------- */
/* File: example2_str.c */
/* ------------------------------------------------------------- */
/* String API. */
/* ------------------------------------------------------------- */
/* Include UNURAN header file. */
#include <unuran.h>
/* ------------------------------------------------------------- */
/* In this example we use a generic distribution object */
/* and sample from this distribution. */
/* */
/* The PDF of our distribution is given by */
/* */
/* / 1 - x*x if |x| <= 1 */
/* f(x) = < */
/* \ 0 otherwise */
/* */
/* We use method TDR (Transformed Density Rejection) which */
/* required a PDF and the derivative of the PDF. */
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double x; /* will hold the random number */
/* Declare UNURAN generator object. */
UNUR_GEN *gen; /* generator object */
/* Create the generator object. */
/* Use a generic continuous distribution. */
/* Choose a method: TDR. */
gen = unur_str2gen(
"distr = cont; pdf=\"1-x*x\"; domain=(-1,1) & method=tdr");
/* It is important to check if the creation of the generator */
/* object was successful. Otherwise `gen' is the NULL pointer */
/* and would cause a segmentation fault if used for sampling. */
if (gen == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* Now you can use the generator object `gen' to sample from */
/* the distribution. Eg.: */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* When you do not need the generator object any more, you */
/* can destroy it. */
unur_free(gen);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
Each method for generating random numbers allows several parameters to be modified. If you do not want to use default values, it is possible to change them. The following example illustrates how to change parameters. For details see Methods.
/* ------------------------------------------------------------- */
/* File: example3.c */
/* ------------------------------------------------------------- */
/* Include UNURAN header file. */
#include <unuran.h>
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double x; /* will hold the random number */
double fparams[2]; /* array for parameters for distribution */
/* Declare the three UNURAN objects. */
UNUR_DISTR *distr; /* distribution object */
UNUR_PAR *par; /* parameter object */
UNUR_GEN *gen; /* generator object */
/* Use a predefined standard distribution: */
/* Gaussian with mean 2. and standard deviation 0.5. */
fparams[0] = 2.;
fparams[1] = 0.5;
distr = unur_distr_normal( fparams, 2 );
/* Choose a method: TDR. */
par = unur_tdr_new(distr);
/* Change some of the default parameters. */
/* We want to use T(x)=log(x) for the transformation. */
unur_tdr_set_c( par, 0. );
/* We want to have the variant with immediate acceptance. */
unur_tdr_set_variant_ia( par );
/* We want to use 10 construction points for the setup */
unur_tdr_set_cpoints ( par, 10, NULL );
/* Create the generator object. */
gen = unur_init(par);
/* Notice that this call has also destroyed the parameter */
/* object `par' as a side effect. */
/* It is important to check if the creation of the generator */
/* object was successful. Otherwise `gen' is the NULL pointer */
/* and would cause a segmentation fault if used for sampling. */
if (gen == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* It is possible to reuse the distribution object to create */
/* another generator object. If you do not need it any more, */
/* it should be destroyed to free memory. */
unur_distr_free(distr);
/* Now you can use the generator object `gen' to sample from */
/* the distribution. Eg.: */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* It is possible with method TDR to truncate the distribution */
/* for an existing generator object ... */
unur_tdr_chg_truncated( gen, -1., 0. );
/* ... and sample again. */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* When you do not need the generator object any more, you */
/* can destroy it. */
unur_free(gen);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
Each method for generating random numbers allows several parameters to be modified. If you do not want to use default values, it is possible to change them. The following example illustrates how to change parameters. For details see Methods.
/* ------------------------------------------------------------- */
/* File: example3_str.c */
/* ------------------------------------------------------------- */
/* String API. */
/* ------------------------------------------------------------- */
/* Include UNURAN header file. */
#include <unuran.h>
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double x; /* will hold the random number */
/* Declare UNURAN generator object. */
UNUR_GEN *gen; /* generator object */
/* Create the generator object. */
/* Use a predefined standard distribution: */
/* Gaussian with mean 2. and standard deviation 0.5. */
/* Choose a method: TDR with parameters */
/* c = 0: use T(x)=log(x) for the transformation; */
/* variant "immediate acceptance"; */
/* number of construction points = 10. */
gen = unur_str2gen(
"normal(2,0.5) & method=tdr; c=0.; variant_ia; cpoints=10");
/* It is important to check if the creation of the generator */
/* object was successful. Otherwise `gen' is the NULL pointer */
/* and would cause a segmentation fault if used for sampling. */
if (gen == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* Now you can use the generator object `gen' to sample from */
/* the distribution. Eg.: */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* It is possible with method TDR to truncate the distribution */
/* for an existing generator object ... */
unur_tdr_chg_truncated( gen, -1., 0. );
/* ... and sample again. */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* When you do not need the generator object any more, you */
/* can destroy it. */
unur_free(gen);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
All generator object use the same default uniform random number generator by default. This can be changed to any generator of your choice such that each generator object has its own random number generator or can share it with some other objects. It is also possible to change the default generator at any time. See Using uniform random number generators, for details.
The following example shows how the uniform random number generator can be set or changed for a generator object. It requires the RNGSTREAMS library to be installed and used. Otherwise the example must be modified accordingly.
/* ------------------------------------------------------------- */
/* File: example_rngstreams.c */
/* ------------------------------------------------------------- */
#ifdef UNURAN_SUPPORTS_RNGSTREAM
/* ------------------------------------------------------------- */
/* This example makes use of the RNGSTREAM library for */
/* for generating uniform random numbers. */
/* (see http://statmath.wu.ac.at/software/RngStreams/) */
/* To compile this example you must have set */
/* ./configure --with-urng-rngstream */
/* (Of course the executable has to be linked against the */
/* RNGSTREAM library.) */
/* ------------------------------------------------------------- */
/* Include UNURAN header files. */
#include <unuran.h>
#include <unuran_urng_rngstreams.h>
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double x; /* will hold the random number */
double fparams[2]; /* array for parameters for distribution */
/* Declare the three UNURAN objects. */
UNUR_DISTR *distr; /* distribution object */
UNUR_PAR *par; /* parameter object */
UNUR_GEN *gen; /* generator object */
/* Declare objects for uniform random number generators. */
UNUR_URNG *urng1, *urng2; /* uniform generator objects */
/* The RNGSTREAMS library sets a package seed. */
unsigned long seed[] = {111u, 222u, 333u, 444u, 555u, 666u};
RngStream_SetPackageSeed(seed);
/* RngStreams only: */
/* Make a object for uniform random number generator. */
/* For details see */
/* http://statmath.wu.ac.at/software/RngStreams/ */
urng1 = unur_urng_rngstream_new("urng-1");
if (urng1 == NULL) exit (EXIT_FAILURE);
/* Use a predefined standard distribution: */
/* Beta with parameters 2 and 3. */
fparams[0] = 2.;
fparams[1] = 3.;
distr = unur_distr_beta( fparams, 2 );
/* Choose a method: TDR. */
par = unur_tdr_new(distr);
/* Set uniform generator in parameter object */
unur_set_urng( par, urng1 );
/* Create the generator object. */
gen = unur_init(par);
/* Notice that this call has also destroyed the parameter */
/* object `par' as a side effect. */
/* It is important to check if the creation of the generator */
/* object was successful. Otherwise `gen' is the NULL pointer */
/* and would cause a segmentation fault if used for sampling. */
if (gen == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* It is possible to reuse the distribution object to create */
/* another generator object. If you do not need it any more, */
/* it should be destroyed to free memory. */
unur_distr_free(distr);
/* Now you can use the generator object `gen' to sample from */
/* the distribution. Eg.: */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* Now we want to switch to a different (independent) stream */
/* of uniform random numbers. */
urng2 = unur_urng_rngstream_new("urng-2");
if (urng2 == NULL) exit (EXIT_FAILURE);
unur_chg_urng( gen, urng2 );
/* ... and sample again. */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* When you do not need the generator object any more, you */
/* can destroy it. */
unur_free(gen);
/* We also should destroy the uniform random number generators.*/
unur_urng_free(urng1);
unur_urng_free(urng2);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
#else
#include <stdio.h>
#include <stdlib.h>
int main(void) {
printf("You must enable the RNGSTREAM library to run this example!\n\n");
exit (77); /* exit code for automake check routines */
}
#endif
/* ------------------------------------------------------------- */
Using Method TDR it is easy to sample pairs of antithetic random variates.
/* ------------------------------------------------------------- */
/* File: example_anti.c */
/* ------------------------------------------------------------- */
#ifdef UNURAN_SUPPORTS_PRNG
/* ------------------------------------------------------------- */
/* This example makes use of the PRNG library for generating */
/* uniform random numbers. */
/* (see http://statmath.wu.ac.at/prng/) */
/* To compile this example you must have set */
/* ./configure --with-urng-prng */
/* (Of course the executable has to be linked against the */
/* PRNG library.) */
/* ------------------------------------------------------------- */
/* Example how to sample from two streams of antithetic random */
/* variates from Gaussian N(2,5) and Gamma(4) distribution, resp.*/
/* ------------------------------------------------------------- */
/* Include UNURAN header files. */
#include <unuran.h>
#include <unuran_urng_prng.h>
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double xn, xg; /* will hold the random number */
double fparams[2]; /* array for parameters for distribution */
/* Declare the three UNURAN objects. */
UNUR_DISTR *distr; /* distribution object */
UNUR_PAR *par; /* parameter object */
UNUR_GEN *gen_normal, *gen_gamma;
/* generator objects */
/* Declare objects for uniform random number generators. */
UNUR_URNG *urng1, *urng2; /* uniform generator objects */
/* PRNG only: */
/* Make a object for uniform random number generator. */
/* For details see http://statmath.wu.ac.at/prng/. */
/* The first generator: Gaussian N(2,5) */
/* uniform generator: We use the Mersenne Twister. */
urng1 = unur_urng_prng_new("mt19937(1237)");
if (urng1 == NULL) exit (EXIT_FAILURE);
/* UNURAN generator object for N(2,5) */
fparams[0] = 2.;
fparams[1] = 5.;
distr = unur_distr_normal( fparams, 2 );
/* Choose method TDR with variant PS. */
par = unur_tdr_new( distr );
unur_tdr_set_variant_ps( par );
/* Set uniform generator in parameter object. */
unur_set_urng( par, urng1 );
/* Set auxilliary uniform random number generator. */
/* We use the default generator. */
unur_use_urng_aux_default( par );
/* Alternatively you can create and use your own auxilliary */
/* uniform random number generator: */
/* UNUR_URNG *urng_aux; */
/* urng_aux = unur_urng_prng_new("tt800"); */
/* if (urng_aux == NULL) exit (EXIT_FAILURE); */
/* unur_set_urng_aux( par, urng_aux ); */
/* Create the generator object. */
gen_normal = unur_init(par);
if (gen_normal == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* Destroy distribution object (gen_normal has its own copy). */
unur_distr_free(distr);
/* The second generator: Gamma(4) with antithetic variates. */
/* uniform generator: We use the Mersenne Twister. */
urng2 = unur_urng_prng_new("anti(mt19937(1237))");
if (urng2 == NULL) exit (EXIT_FAILURE);
/* UNURAN generator object for gamma(4) */
fparams[0] = 4.;
distr = unur_distr_gamma( fparams, 1 );
/* Choose method TDR with variant PS. */
par = unur_tdr_new( distr );
unur_tdr_set_variant_ps( par );
/* Set uniform generator in parameter object. */
unur_set_urng( par, urng2 );
/* Set auxilliary uniform random number generator. */
/* We use the default generator. */
unur_use_urng_aux_default( par );
/* Alternatively you can create and use your own auxilliary */
/* uniform random number generator (see above). */
/* Notice that both generator objects gen_normal and */
/* gen_gamma can share the same auxilliary URNG. */
/* Create the generator object. */
gen_gamma = unur_init(par);
if (gen_gamma == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* Destroy distribution object (gen_normal has its own copy). */
unur_distr_free(distr);
/* Now we can sample pairs of negatively correlated random */
/* variates. E.g.: */
for (i=0; i<10; i++) {
xn = unur_sample_cont(gen_normal);
xg = unur_sample_cont(gen_gamma);
printf("%g, %g\n",xn,xg);
}
/* When you do not need the generator objects any more, you */
/* can destroy it. */
unur_free(gen_normal);
unur_free(gen_gamma);
/* We also should destroy the uniform random number generators.*/
unur_urng_free(urng1);
unur_urng_free(urng2);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
#else
#include <stdio.h>
#include <stdlib.h>
int main(void) {
printf("You must enable the PRNG library to run this example!\n\n");
exit (77); /* exit code for automake check routines */
}
#endif
/* ------------------------------------------------------------- */
Using Method TDR it is easy to sample pairs of antithetic random variates.
/* ------------------------------------------------------------- */
/* File: example_anti_str.c */
/* ------------------------------------------------------------- */
/* String API. */
/* ------------------------------------------------------------- */
#ifdef UNURAN_SUPPORTS_PRNG
/* ------------------------------------------------------------- */
/* This example makes use of the PRNG library for generating */
/* uniform random numbers. */
/* (see http://statmath.wu.ac.at/prng/) */
/* To compile this example you must have set */
/* ./configure --with-urng-prng */
/* (Of course the executable has to be linked against the */
/* PRNG library.) */
/* ------------------------------------------------------------- */
/* Example how to sample from two streams of antithetic random */
/* variates from Gaussian N(2,5) and Gamma(4) distribution, resp.*/
/* ------------------------------------------------------------- */
/* Include UNURAN header files. */
#include <unuran.h>
#include <unuran_urng_prng.h>
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double xn, xg; /* will hold the random number */
/* Declare UNURAN generator objects. */
UNUR_GEN *gen_normal, *gen_gamma;
/* PRNG only: */
/* Make a object for uniform random number generator. */
/* For details see http://statmath.wu.ac.at/prng/. */
/* Create the first generator: Gaussian N(2,5) */
gen_normal = unur_str2gen("normal(2,5) & method=tdr; variant_ps & urng=mt19937(1237)");
if (gen_normal == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* Set auxilliary uniform random number generator. */
/* We use the default generator. */
unur_chgto_urng_aux_default(gen_normal);
/* The second generator: Gamma(4) with antithetic variates. */
gen_gamma = unur_str2gen("gamma(4) & method=tdr; variant_ps & urng=anti(mt19937(1237))");
if (gen_gamma == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
unur_chgto_urng_aux_default(gen_gamma);
/* Now we can sample pairs of negatively correlated random */
/* variates. E.g.: */
for (i=0; i<10; i++) {
xn = unur_sample_cont(gen_normal);
xg = unur_sample_cont(gen_gamma);
printf("%g, %g\n",xn,xg);
}
/* When you do not need the generator objects any more, you */
/* can destroy it. */
/* But first we have to destroy the uniform random number */
/* generators. */
unur_urng_free(unur_get_urng(gen_normal));
unur_urng_free(unur_get_urng(gen_gamma));
unur_free(gen_normal);
unur_free(gen_gamma);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
#else
#include <stdio.h>
#include <stdlib.h>
int main(void) {
printf("You must enable the PRNG library to run this example!\n\n");
exit (77); /* exit code for automake check routines */
}
#endif
/* ------------------------------------------------------------- */
See Methods for continuous univariate distributions.
See Methods for continuous empirical univariate distributions.
See Methods for continuous empirical multivariate distributions.
See Methods for discrete univariate distributions.
The string interface (string API) provided by the
unur_str2gen
call is the easiest way to use UNU.RAN. This
function takes a character string as its argument. The string is
parsed and the information obtained is used to create a generator
object. It returns NULL if this fails, either due to a syntax
error, or due to invalid data. In both cases unur_error is
set to the corresponding error codes
(see Error reporting).
Additionally there exists the call
unur_str2distr
that only
produces a distribution object.
Notice that the string interface does not implement all features of the UNU.RAN library. For trickier tasks it might be necessary to use the UNU.RAN calls.
In Examples, all examples are given using both the UNU.RAN standard API and this convenient string API. The corresponding programm codes are equivalent.
Get a generator object for the distribution, method and uniform random number generator as described in the given string. See Syntax of String Interface, for details.
Get a distribution object for the distribution described in string. See Syntax of String Interface, and Distribution String, for details. However, only the block for the distribution object is allowed.
— : UNUR_GEN* unur_makegen_dsu (const UNUR_DISTR* distribution, const char* methodstr, UNUR_URNG* urng)
Make a generator object for the distribution, method and uniform random number generator. The distribution can be given either as string distrstr or as a distribution object distr. The method must be given as a string methodstr. For the syntax of these strings see Syntax of String Interface. However, the
methodkeyword is optional for these calls and can be omitted. If methodstr is the empty (blank) string orNULLmethod AUTO is used. The uniform random number generator is optional. If urng isNULLthen the default uniform random number generator is used.
The given string holds information about the requested distribution and (optional) about the sampling method and the uniform random number generator invoked. The interpretation of the string is not case-sensitive, all white spaces are ignored.
The string consists of up to three blocks, separated by ampersands
&.
Each block consists of <key>=<value> pairs, separated by
semicolons ;.
The first key in each block is used to indicate each block. We have three different blocks with the following (first) keys:
distrmethodurngThe distr block must be the very first block and is
obligatory. All the other blocks are optional and can be arranged
in arbitrary order.
For details see the following description of each block.
In the following example
distr = normal(3.,0.75); domain = (0,inf) & method = tdr; c = 0
we have a distribution block for the truncated normal distribution with mean 3 and standard deviation 0.75 on domain (0,infinity); and block for choosing method TDR with parameter c set to 0.
The
<key>=<value> pairs that follow the first (initial) pair
in each block are used to set parameters.
The name of the parameter is given by the <key> string. It is
deduced from the UNU.RAN set calls by taking the part after
..._set_.
The <value> string holds the parameters to be
set, separated by commata ,.
There are three types of parameters:
"...""...",
(...,...),,
enclosed in parenthesis (...), and
"..." or parenthesis (...).
It is interpreted as float or integer depending on the type of
the corresponding parameter.
<value> string (including the character =) can be
omitted when no argument is required.
At the moment not all set calls are supported.
The syntax for the <value> can be directly derived from the
corresponding set calls. To simplify the syntax additional
shortcuts are possible. The following table lists the parameters for
the set calls that are supported by the string interface; the
entry in parenthesis gives the type of the argument as
<value> string:
int (number):true and on are transformed to 1,
false and off are transformed to 0.
A missing argument is interpreted as 1.
int, int (number, number or list):inf and -inf are transformed to INT_MAX and
INT_MIN respectively, i.e. the largest and smallest
integers that can be represented by the computer.
unsigned (number):double (number):inf is transformed to UNUR_INFINITY.
double, double (number, number or list):inf is transformed to UNUR_INFINITY. However using
inf in the list might not work for all versions of C. Then it
is recommended to use two single numbers instead of a list.
int, double* ([number,] list or number):NULL
pointer is used instead an array as argument.
double*, int (list [,number]):distribution block!)
char* (string):Notice that missing entries in a list of numbers are interpreted as
0. E.g, a the list (1,,3) is read as (1,0,3), the
list (1,2,) as (1,2,0).
The the list of key strings in
Keys for Distribution String, and
Keys for Method String, for further details.
The distr block must be the very first block and is
obligatory. For that reason the keyword distr is optional and
can be omitted (together with the = character).
Moreover it is ignored while parsing the string. However, to
avoid some possible confusion it has to start with the
letter d (if it is given at all).
The value of the distr key is used to get the distribution
object, either via a unur_distr_<value> call for a standard
distribution via a unur_distr_<value>_new call to get an
object of a generic distribution.
However not all generic distributions are supported yet.
The parameters for the standard distribution are given
as a list. There must not be any character (other than white space)
between the name of the standard distribution and the opening
parenthesis ( of this list. E.g., to get a beta distribution,
use
distr = beta(2,4)
To get an object for a discrete distribution with probability vector (0.5,0.2,0.3), use
distr = discr; pv = (0.5,0.2,0.3)
It is also possible to set a PDF, PMF, or CDF using a string.
E.g., to create a continuous distribution with PDF proportional to
exp(-sqrt(2+(x-1)^2) + (x-1)) and domain (0,inf) use
distr = cont; pdf = "exp(-sqrt(2+(x-1)^2) + (x-1))"
Notice: If this string is used in an
unur_str2distr
or
unur_str2gen
call the double quotes " must be protected by
\". Alternatively, single quotes may be used instead
distr = cont; pdf = 'exp(-sqrt(2+(x-1)^2) + (x-1))'
For the details of function strings see Function String.
List of standard distributions see Standard distributions
[distr =] beta(...) ⇒ see beta
[distr =] binomial(...) ⇒ see binomial
[distr =] cauchy(...) ⇒ see cauchy
[distr =] chi(...) ⇒ see chi
[distr =] chisquare(...) ⇒ see chisquare
[distr =] exponential(...) ⇒ see exponential
[distr =] extremeI(...) ⇒ see extremeI
[distr =] extremeII(...) ⇒ see extremeII
[distr =] F(...) ⇒ see F
[distr =] gamma(...) ⇒ see gamma
[distr =] geometric(...) ⇒ see geometric
[distr =] gig(...) ⇒ see gig
[distr =] gig2(...) ⇒ see gig2
[distr =] hyperbolic(...) ⇒ see hyperbolic
[distr =] hypergeometric(...) ⇒ see hypergeometric
[distr =] ig(...) ⇒ see ig
[distr =] laplace(...) ⇒ see laplace
[distr =] logarithmic(...) ⇒ see logarithmic
[distr =] logistic(...) ⇒ see logistic
[distr =] lognormal(...) ⇒ see lognormal
[distr =] lomax(...) ⇒ see lomax
[distr =] negativebinomial(...) ⇒ see negativebinomial
[distr =] normal(...) ⇒ see normal
[distr =] pareto(...) ⇒ see pareto
[distr =] poisson(...) ⇒ see poisson
[distr =] powerexponential(...) ⇒ see powerexponential
[distr =] rayleigh(...) ⇒ see rayleigh
[distr =] slash(...) ⇒ see slash
[distr =] student(...) ⇒ see student
[distr =] triangular(...) ⇒ see triangular
[distr =] uniform(...) ⇒ see uniform
[distr =] weibull(...) ⇒ see weibull
List of generic distributions see Handling Distribution Objects
Notice: Order statistics for continuous distributions (see CORDER) are supported by using the key
orderstatistics for distributions
of type CONT.
List of keys that are available via the String API. For description see the corresponding UNU.RAN set calls.
name = "<string>"
cemp (Distribution Type) (see CEMP)
data = (<list>) [, <int>]hist_bins = (<list>) [, <int>]hist_domain = <double>, <double> | (<list>)hist_prob = (<list>) [, <int>]
cont (Distribution Type) (see CONT)
cdf = "<string>"center = <double>domain = <double>, <double> | (<list>)hr = "<string>"logcdf = "<string>"logpdf = "<string>"mode = <double>pdf = "<string>"pdfarea = <double>pdfparams = (<list>) [, <int>]orderstatistics = <int>, <int> | (<list>)
discr (Distribution Type) (see DISCR)
cdf = "<string>"domain = <int>, <int> | (<list>)mode [= <int>]pmf = "<string>"pmfparams = (<list>) [, <int>]pmfsum = <double>pv = (<list>) [, <int>]
In unuran it is also possible to define functions (e.g. CDF or PDF) as strings. As you can see in Example 2 (Example_2_str) it is very easy to define the PDF of a distribution object by means of a string. The possibilities using this string interface are more restricted than using a pointer to a routine coded in C (Example_2). But the differences in evaluation time is small. When a distribution object is defined using this string interface then of course the same conditions on the given density or CDF must be satisfied for a chosen method as for the standard API. This string interface can be used for both within the UNU.RAN string API using the unur_str2gen call, and for calls that define the density or CDF for a particular distribution object as done with (e.g.) the call unur_distr_cont_set_pdfstr. Here is an example for the latter case:
unur_distr_cont_set_pdfstr(distr,"1-x*x");
The syntax for the function string is case insensitive, white spaces are ingnored. The expressions are similar to most programming languages and mathematical programs (see also the examples below). It is especially influenced by C. The usual preceedence rules are used (from highest to lowest preceedence: functions, power, multiplication, addition, relation operators). Use parentheses in case of doubt or when these preceedences should be changed.
Relation operators can be used as indicator functions, i.e. the term
(x>1) is evaluted as 1 if this relation is satisfied,
and as 0 otherwise.
The first unknown symbol (letter or word) is interpreted as the
variable of the function. It is recommended to use x.
Only one variable can be used.
Important: The symbol
e is used twice, for Euler's
constant (= 2.7182...) and as exponent.
The multiplication operator * must not be omitted, i.e.
2 x is interpreted as the string 2x (which will result
in a syntax error).
Numbers
Numbers are composed using digits and, optionally, a sign, a decimal point, and an exponent indicated by e.
|
Constants
|
Special symbols
|
Relation operators (Indicator functions)
|
Arithmetic operators
|
Functions
|
Variable
|
1.231+7.9876*x-1.234e-3*x^2+3.335e-5*x^3
sin(2*pi*x)+x^2
exp(-((x-3)/2.1)^2)
It is also possible to define functions using different terms
on separate domains. However, instead of constructs using
if ... then ... else ...
indicator functions are available.
For example to define the density of triangular distribution with domain (-1,1) and mode 0 use
(x>-1)*(x<0)*(1+x) + (x>=0)*(x<1)*(1-x)
The key method is obligatory, it must be the first key and its
value is the name of a method suitable for the choosen standard
distribution. E.g., if method AROU is chosen, use
method = arou
Of course the all following keys dependend on the method choosen at
first. All corresponding set calls of UNU.RAN are available
and the key is the string after the unur_<methodname>_set_
part of the command. E.g., UNU.RAN provides the command
unur_arou_set_max_sqhratio to set a parameter of method AROU.
To call this function via the string-interface, the
key max_sqhratio can be used:
max_sqhratio = 0.9
Additionally the keyword debug can be used to set debugging
flags (see Debugging, for details).
If this block is omitted, a suitable default method is used. Notice however that the default method may change in future versions of UNU.RAN.
List of methods and keys that are available via the String API. For description see the corresponding UNU.RAN set calls.
method = arou ⇒ unur_arou_new
(see AROU)
cpoints = <int> [, (<list>)] | (<list>)darsfactor = <double>guidefactor = <double>max_segments [= <int>]max_sqhratio = <double>pedantic [= <int>]usecenter [= <int>]usedars [= <int>]verify [= <int>]
method = ars ⇒ unur_ars_new
(see ARS)
cpoints = <int> [, (<list>)] | (<list>)max_intervals [= <int>]max_iter [= <int>]pedantic [= <int>]reinit_ncpoints [= <int>]reinit_percentiles = <int> [, (<list>)] | (<list>)verify [= <int>]
method = auto ⇒ unur_auto_new
(see AUTO)
logss [= <int>]
method = cstd ⇒ unur_cstd_new
(see CSTD)
variant = <unsigned>
method = dari ⇒ unur_dari_new
(see DARI)
cpfactor = <double>squeeze [= <int>]tablesize [= <int>]verify [= <int>]
method = dau ⇒ unur_dau_new
(see DAU)
urnfactor = <double>
method = dgt ⇒ unur_dgt_new
(see DGT)
guidefactor = <double>variant = <unsigned>
method = dsrou ⇒ unur_dsrou_new
(see DSROU)
cdfatmode = <double>verify [= <int>]
method = dstd ⇒ unur_dstd_new
(see DSTD)
variant = <unsigned>
method = empk ⇒ unur_empk_new
(see EMPK)
beta = <double>kernel = <unsigned>positive [= <int>]smoothing = <double>varcor [= <int>]
method = gibbs ⇒ unur_gibbs_new
(see GIBBS)
burnin [= <int>]c = <double>thinning [= <int>]variant_coordinatevariant_random_direction
method = hinv ⇒ unur_hinv_new
(see HINV)
boundary = <double>, <double> | (<list>)cpoints = (<list>), <int>guidefactor = <double>max_intervals [= <int>]order [= <int>]u_resolution = <double>
method = hitro ⇒ unur_hitro_new
(see HITRO)
adaptive_multiplier = <double>burnin [= <int>]r = <double>thinning [= <int>]use_adaptiveline [= <int>]use_adaptiverectangle [= <int>]use_boundingrectangle [= <int>]v = <double>variant_coordinatevariant_random_direction
method = hrb ⇒ unur_hrb_new
(see HRB)
upperbound = <double>verify [= <int>]
method = hrd ⇒ unur_hrd_new
(see HRD)
verify [= <int>]
method = hri ⇒ unur_hri_new
(see HRI)
p0 = <double>verify [= <int>]
method = itdr ⇒ unur_itdr_new
(see ITDR)
cp = <double>ct = <double>verify [= <int>]xi = <double>
method = mvtdr ⇒ unur_mvtdr_new
(see MVTDR)
boundsplitting = <double>maxcones [= <int>]stepsmin [= <int>]verify [= <int>]
method = ninv ⇒ unur_ninv_new
(see NINV)
max_iter [= <int>]start = <double>, <double> | (<list>)table [= <int>]u_resolution = <double>usebisectusenewtonuseregulax_resolution = <double>
method = nrou ⇒ unur_nrou_new
(see NROU)
center = <double>r = <double>u = <double>, <double> | (<list>)v = <double>verify [= <int>]
method = pinv ⇒ unur_pinv_new
(see PINV)
boundary = <double>, <double> | (<list>)max_intervals [= <int>]order [= <int>]searchboundary = <int>, <int> | (<list>)u_resolution = <double>usecdfusepdf
method = srou ⇒ unur_srou_new
(see SROU)
cdfatmode = <double>pdfatmode = <double>r = <double>usemirror [= <int>]usesqueeze [= <int>]verify [= <int>]
method = ssr ⇒ unur_ssr_new
(see SSR)
cdfatmode = <double>pdfatmode = <double>usesqueeze [= <int>]verify [= <int>]
method = tabl ⇒ unur_tabl_new
(see TABL)
areafraction = <double>boundary = <double>, <double> | (<list>)cpoints = <int> [, (<list>)] | (<list>)darsfactor = <double>guidefactor = <double>max_intervals [= <int>]max_sqhratio = <double>nstp [= <int>]pedantic [= <int>]slopes = (<list>), <int>usedars [= <int>]useear [= <int>]variant_ia [= <int>]variant_splitmode = <unsigned>verify [= <int>]
method = tdr ⇒ unur_tdr_new
(see TDR)
c = <double>cpoints = <int> [, (<list>)] | (<list>)darsfactor = <double>guidefactor = <double>max_intervals [= <int>]max_sqhratio = <double>pedantic [= <int>]reinit_ncpoints [= <int>]reinit_percentiles = <int> [, (<list>)] | (<list>)usecenter [= <int>]usedars [= <int>]usemode [= <int>]variant_gwvariant_iavariant_psverify [= <int>]
method = utdr ⇒ unur_utdr_new
(see UTDR)
cpfactor = <double>deltafactor = <double>pdfatmode = <double>verify [= <int>]
method = vempk ⇒ unur_vempk_new
(see VEMPK)
smoothing = <double>varcor [= <int>]
method = vnrou ⇒ unur_vnrou_new
(see VNROU)
r = <double>v = <double>verify [= <int>]
The value of the urng key is passed to the PRNG interface (see
PRNG manual
for details).
However it only works when using the PRNG library is enabled,
see Installation for details. There are no other keys.
IMPORTANT: UNU.RAN creates a new uniform random number generator for the generator object. The pointer to this uniform generator has to be read and saved via a unur_get_urng call in order to clear the memory before the UNU.RAN generator object is destroyed.
If this block is omitted the UNU.RAN default generator is used (which must not be destroyed).
Objects of type UNUR_DISTR are used for handling
distributions. All data about a distribution are stored in this
object. UNU.RAN provides functions that return instances of such
objects for standard distributions
(see Standard distributions).
It is then possible to change these distribution objects by
various set calls. Moreover, it is possible to build a
distribution object entirely from scratch. For this purpose
there exists unur_distr_<type>_new calls that
return an empty object of this type for each object type
(eg. univariate contiuous) which can be filled with the
appropriate set calls.
UNU.RAN distinguishes between several types of distributions, each of which has its own sets of possible parameters (for details see the corresponding sections):
Notice that there are essential data about a distribution, eg. the PDF, a list of (shape, scale, location) parameters for the distribution, and the domain of (the possibly truncated) distribution. And there exist parameters that are/can be derived from these, eg. the mode of the distribution or the area below the given PDF (which need not be normalized for many methods). UNU.RAN keeps track of parameters which are known. Thus if one of the essential parameters is changed all derived parameters are marked as unknown and must be set again if these are required for the chosen generation method. Additionally to set calls there are calls for updating derived parameters for objects provided by the UNU.RAN library of standard distributions (one for each parameter to avoid computational overhead since not all parameters are required for all generator methods).
All parameters of distribution objects can be read by corresponding get calls.
Every generator object has its own copy of a distribution object which is accessible by a unur_get_distr call. Thus the parameter for this distribution can be read. However, never extract the distribution object out of a generator object and run one of the set calls on it to modify the distribution. (How should the poor generator object know what has happend?) Instead there exist calls for each of the generator methods that change particular parameters of the internal copy of the distribution object.
UNU.RAN collects all data required for a particular generation method in a distribution object. There are two ways to get an instance of a distributions object:
In either cases the corresponding
unur_distr_<type>_set_<param> calls to set the
necessary parameters <param> (case 1), or
change the values of the standard distribution in case 2 (if
this makes sense for you). In the latter case <type>
is the type to which the standard distribution belongs to.
These set calls return UNUR_SUCCESS when the
correspondig parameter has been set successfully. Otherwise an
error code is returned.
The parameters of a distribution are divided into essential and derived parameters.
Notice, that there are some restrictions in setting parameters
to avoid possible confusions.
Changing essential parameters marks derived parameters as
unknown. Some of the parameters cannot be changed any
more when already set; some parameters block each others.
In such a case a new instance of a distribution object has to be
build.
Additionally unur_distr_<type>_upd_<param> calls can be used for updating derived parameters for objects provided by the UNU.RAN library of standard distributions.
All parameters of a distribution object get be read by means of unur_distr_<type>_get_<param> calls.
Every distribution object be identified by its name which
is a string of arbitrary characters provided by the user. For
standard distribution it is automatically set to
<name> in the corresponding new call. It can
be changed to any other string.
The calls in this section can be applied to all distribution objects.
— : const char* unur_distr_get_name (const UNUR_DISTR* distribution)
Set and get name of distribution. The name can be an arbitrary character string. It can be used to identify generator objects for the user. It is used by UNU.RAN when printing information of the distribution object into a log files.
Get number of components of a random vector (its dimension) the distribution.
For univariate distributions it returns dimension
1.For matrix distributions it returns the number of components (i.e., number of rows times number of columns). When the respective numbers of rows and columns are needed use unur_distr_matr_get_dim instead.
Get type of distribution. Possible types are
UNUR_DISTR_CONT- univariate continuous distribution
UNUR_DISTR_CEMP- empirical continuous univariate distribution (i.e. a sample)
UNUR_DISTR_CVEC- continuous mulitvariate distribution
UNUR_DISTR_CVEMP- empirical continuous multivariate distribution (i.e. a vector sample)
UNUR_DISTR_DISCR- discrete univariate distribution
UNUR_DISTR_MATR- matrix distribution
Alternatively the unur_distr_is_<TYPE> calls can be used.
TRUEif distribution is a continuous univariate distribution.
TRUEif distribution is a continuous multivariate distribution.
TRUEif distribution is an empirical continuous univariate distribution, i.e. a sample.
TRUEif distribution is an empirical continuous multivariate distribution.
TRUEif distribution is a discrete univariate distribution.
TRUEif distribution is a matrix distribution.
Store a pointer to an external object. This might be usefull if the PDF, PMF, CDF or other functions used to implement a particular distribution a parameter set that cannot be stored as doubles (e.g. pointers to some structure that holds information of the distribution).
Important: When UNU.RAN copies this distribution object into the generator object, then the address extobj that this pointer contains is simply copied. Thus the generator holds an address of a non-private object! Once the generator object has been created any change in the external object might effect the generator object.
Warning: External objects must be used with care. Once the generator object has been created or the distribution object has been copied you must not destroy this external object.
Get the pointer to the external object.
Important: Changing this object must be done with with extreme care.
The calls in this section can be applied to continuous univariate distributions.
Some methods use the logarithm of the density if available.
Create a new (empty) object for univariate continuous distribution.
— : int unur_distr_cont_set_dpdf (UNUR_DISTR* distribution, UNUR_FUNCT_CONT* dpdf)
— : int unur_distr_cont_set_cdf (UNUR_DISTR* distribution, UNUR_FUNCT_CONT* cdf)
— : int unur_distr_cont_set_invcdf (UNUR_DISTR* distribution, UNUR_FUNCT_CONT* invcdf)
Set respective pointer to the probability density function (PDF), the derivative of the probability density function (dPDF), the cumulative distribution function (CDF), and the inverse CDF of the distribution. Each of these function pointers must be of type
double funct(double x, const UNUR_DISTR *distr).Due to the fact that some of the methods do not require a normalized PDF the following is important:
- The given CDF must be the cumulative distribution function of the (non-truncated) distribution. If a distribution from the UNU.RAN library of standard distributions (see Standard distributions) is truncated, there is no need to change the CDF.
- If both the CDF and the PDF are used (for a method or for order statistics), the PDF must be the derivative of the CDF. If a truncated distribution for one of the standard distributions from the UNU.RAN library of standard distributions is used, there is no need to change the PDF.
- If the area below the PDF is required for a given distribution it must be given by the unur_distr_cont_set_pdfarea call. For a truncated distribution this must be of course the integral of the PDF in the given truncated domain. For distributions from the UNU.RAN library of standard distributions this is done automatically by the unur_distr_cont_upd_pdfarea call.
It is important to note that all these functions must return a result for all values of x. Eg., if the domain of a given PDF is the interval [-1,1], then the given function must return
0.0for all points outside this interval. In case of an overflow the PDF should returnUNUR_INFINITY.It is not possible to change such a function. Once the PDF or CDF is set it cannot be overwritten. This also holds when the logPDF is given or when the PDF is given by the unur_distr_cont_set_pdfstr or unur_distr_cont_set_logpdfstr call. A new distribution object has to be used instead.
— : UNUR_FUNCT_CONT* unur_distr_cont_get_dpdf (const UNUR_DISTR* distribution)
— : UNUR_FUNCT_CONT* unur_distr_cont_get_cdf (const UNUR_DISTR* distribution)
— : UNUR_FUNCT_CONT* unur_distr_cont_get_invcdf (const UNUR_DISTR* distribution)
Get the respective pointer to the PDF, the derivative of the PDF, the CDF, and the inverse CDF of the distribution. The pointer is of type
double funct(double x, const UNUR_DISTR *distr). If the corresponding function is not available for the distribution, theNULLpointer is returned.
— : double unur_distr_cont_eval_dpdf (double x, const UNUR_DISTR* distribution)
— : double unur_distr_cont_eval_cdf (double x, const UNUR_DISTR* distribution)
— : double unur_distr_cont_eval_invcdf (double u, const UNUR_DISTR* distribution)
Evaluate the PDF, derivative of the PDF, the CDF, and the inverse CDF at x and u,respectively. Notice that distribution must not be the
NULLpointer. If the corresponding function is not available for the distribution,UNUR_INFINITYis returned andunur_errnois set toUNUR_ERR_DISTR_DATA.IMPORTANT: In the case of a truncated standard distribution these calls always return the respective values of the untruncated distribution!
— : int unur_distr_cont_set_dlogpdf (UNUR_DISTR* distribution, UNUR_FUNCT_CONT* dlogpdf)
— : int unur_distr_cont_set_logcdf (UNUR_DISTR* distribution, UNUR_FUNCT_CONT* logcdf)
— : UNUR_FUNCT_CONT* unur_distr_cont_get_logpdf (const UNUR_DISTR* distribution)
— : UNUR_FUNCT_CONT* unur_distr_cont_get_dlogpdf (const UNUR_DISTR* distribution)
— : UNUR_FUNCT_CONT* unur_distr_cont_get_logcdf (const UNUR_DISTR* distribution)
— : double unur_distr_cont_eval_logpdf (double x, const UNUR_DISTR* distribution)
— : double unur_distr_cont_eval_dlogpdf (double x, const UNUR_DISTR* distribution)
— : double unur_distr_cont_eval_logcdf (double x, const UNUR_DISTR* distribution)
Analogous calls for the logarithm of the density distribution functions.
This function provides an alternative way to set a PDF and its derivative of the distribution. pdfstr is a character string that contains the formula for the PDF, see Function String, for details. The derivative of the given PDF is computed automatically. See also the remarks for the unur_distr_cont_set_pdf call.
It is not possible to call this funtion twice or to call this function after a unur_distr_cont_set_pdf call.
This function provides an alternative way to set a CDF; analogously to the unur_distr_cont_set_pdfstr call. The PDF and its derivative of the given CDF are computed automatically.
— : char* unur_distr_cont_get_dpdfstr (const UNUR_DISTR* distribution)
— : char* unur_distr_cont_get_cdfstr (const UNUR_DISTR* distribution)
Get pointer to respective string for PDF, derivate of PDF, and CDF of distribution that is given as string (instead of a function pointer). This call allocates memory to produce this string. It should be freed when it is not used any more.
Sets array of parameters for distribution. There is an upper limit for the number of parameters
n_params. It is given by the macroUNUR_DISTR_MAXPARAMSin unuran_config.h. (It is set to 5 by default but can be changed to any appropriate nonnegative number.) If n_params is negative or exceeds this limit no parameters are copied into the distribution object andunur_errnois set toUNUR_ERR_DISTR_NPARAMS.For standard distributions from the UNU.RAN library the parameters are checked. Moreover, the domain is updated automatically unless it has been changed before by a unur_distr_cont_set_domain call. If the given parameters are invalid for the standard distribution, then no parameters are set and an error code is returned. Notice, that the given parameter list for such a distribution is handled in the same way as in the corresponding new calls, i.e. optional parameters for the PDF that are not present in the given list are (re-)set to their default values.
Important: If the parameters of a distribution from the UNU.RAN library of standard distributions (see Standard distributions) are changed, then neither its mode nor the normalization constant are updated. Please use the respective calls unur_distr_cont_upd_mode and unur_distr_cont_upd_pdfarea. Moreover, if the domain has been changed by a unur_distr_cont_set_domain it is not automatically updated, either. Updating the normalization constant is in particular very important, when the CDF of the distribution is used.
Get number of parameters of the PDF and set pointer params to array of parameters. If no parameters are stored in the object, an error code is returned and
paramsis set toNULL.Important: Do not change the entries in params!
This function provides an interface for additional vector parameters for a continuous distribution.
It sets the parameter with number par. par indicates directly which of the parameters is set and must be a number between
0andUNUR_DISTR_MAXPARAMS-1 (the upper limit of possible parameters defined in unuran_config.h; it is set to 5 but can be changed to any appropriate nonnegative number.)The entries of a this parameter are given by the array param_vec of size n_param_vec.
If param_vec is
NULLthen the corresponding entry is cleared.If an error occurs no parameters are copied into the parameter object
unur_errnois set toUNUR_ERR_DISTR_DATA.
Get parameter of the PDF with number par. The pointer to the parameter array is stored in param_vecs, its size is returned by the function. If the requested parameter is not set, then an error code is returned and
paramsis set toNULL.Important: Do not change the entries in param_vecs!
— : char* unur_distr_cont_get_logpdfstr (const UNUR_DISTR* distribution)
— : char* unur_distr_cont_get_dlogpdfstr (const UNUR_DISTR* distribution)
— : int unur_distr_cont_set_logcdfstr (UNUR_DISTR* distribution, const char* logcdfstr)
— : char* unur_distr_cont_get_logcdfstr (const UNUR_DISTR* distribution)
Analogous calls for the logarithm of the density and distribution functions.
Set the left and right borders of the domain of the distribution. This can also be used to truncate an existing distribution. For setting the boundary to +/- infinityuse
+/- UNUR_INFINITY. If right is not strictly greater than left no domain is set andunur_errnois set toUNUR_ERR_DISTR_SET.Important: For some technical reasons it is assumed that the density is unimodal and thus monotone on either side of the mode! This is used in the case when the given mode is outside of the original domain. Then the mode is set to the corresponding boundary of the new domain. If this result is not the desired it must be changed by using a unur_distr_cont_set_mode call (or a unur_distr_cont_upd_mode call). The same holds for the center of the distribution.
Get the left and right borders of the domain of the distribution. If the domain is not set
+/- UNUR_INFINITYis assumed and returned. No error is reported in this case.
Get the left and right borders of the (truncated) domain of the distribution. For non-truncated distribution this call is equivalent to the unur_distr_cont_get_domain call.
This call is only useful in connection with a unur_get_distr call to get the boundaries of the sampling region of a generator object.
Set pointer to the hazard rate (HR) of the distribution.
The hazard rate (or failure rate) is a mathematical way of describing aging. If the lifetime X is a random variable with density f(x) and CDF F(x) the hazard rate h(x) is defined as h(x) = f(x) / (1-F(x)). In other words, h(x) represents the (conditional) rate of failure of a unit that has survived up to time x with probability 1-F(x). The key distribution is the exponential distribution as it has constant hazard rate of value 1. Hazard rates tending to infinity describe distributions with sub-exponential tails whereas distributions with hazard rates tending to zero have heavier tails than the exponential distribution.
It is important to note that all these functions must return a result for all floats x. In case of an overflow the PDF should return
UNUR_INFINITY.Important: Do not simply use f(x) / (1-F(x)), since this is numerically very unstable and results in numerical noise if F(x) is (very) close to 1. Moreover, if the density f(x) is known a generation method that uses the density is more appropriate.
It is not possible to change such a function. Once the HR is set it cannot be overwritten. This also holds when the HR is given by the unur_distr_cont_set_hrstr call. A new distribution object has to be used instead.
Get the pointer to the hazard rate of the distribution. The pointer is of type
double funct(double x, const UNUR_DISTR *distr). If the corresponding function is not available for the distribution, theNULLpointer is returned.
Evaluate the hazard rate at x. Notice that distribution must not be the
NULLpointer. If the corresponding function is not available for the distribution,UNUR_INFINITYis returned andunur_errnois set toUNUR_ERR_DISTR_DATA.
This function provides an alternative way to set a hazard rate and its derivative of the distribution. hrstr is a character string that contains the formula for the HR, see Function String, for details. See also the remarks for the unur_distr_cont_set_hr call.
It is not possible to call this funtion twice or to call this function after a unur_distr_cont_set_hr call.
Get pointer to string for HR of distribution that is given via the string interface. This call allocates memory to produce this string. It should be freed when it is not used any more.
The following paramters must be set whenever one of the essential parameters has been set or changed (and the parameter is required for the chosen method).
Set mode of distribution. The mode must be contained in the domain of distribution. Otherwise the mode is not set and
unur_errnois set toUNUR_ERR_DISTR_SET. For distributions with unbounded density, this call is used to set the pole of the PDF. Notice that the PDF should then return UNUR_INFINITY at the pole. Notice that the mode is adjusted when the domain is set, see the remark for the unur_distr_cont_set_domain call.
Recompute the mode of the distribution. This call works properly for distribution objects from the UNU.RAN library of standard distributions when the corresponding function is available. Otherwise a (slow) numerical mode finder based on Brent's algorithm is used. If it failes
unur_errnois set toUNUR_ERR_DISTR_DATA.
Get mode of distribution. If the mode is not marked as known, unur_distr_cont_upd_mode is called to compute the mode. If this is not successful
UNUR_INFINITYis returned andunur_errnois set toUNUR_ERR_DISTR_GET. (There is no difference between the case where no routine for computing the mode is available and the case where no mode exists for the distribution at all.)
Set center of the distribution. The center is used by some methods to shift the distribution in order to decrease numerical round-off error. If not given explicitly a default is used.
Important: This call does not check whether the center is contained in the given domain.
Default: The mode, if set by a unur_distr_cont_set_mode or unur_distr_cont_upd_mode call; otherwise
0.
Get center of the distribution. It always returns some point as there always exists a default for the center, see unur_distr_cont_set_center.
Set the area below the PDF. If
areais non-positive, no area is set andunur_errnois set toUNUR_ERR_DISTR_SET.For a distribution object created by the UNU.RAN library of standard distributions you always should use the unur_distr_cont_upd_pdfarea. Otherwise there might be ambiguous side-effects.
Recompute the area below the PDF of the distribution. It only works for distribution objects from the UNU.RAN library of standard distributions when the corresponding function is available. Otherwise
unur_errnois set toUNUR_ERR_DISTR_DATA.This call also sets the normalization constant such that the given PDF is the derivative of a given CDF, i.e. the area is 1. However, for truncated distributions the area is smaller than 1.
The call does not work for distributions from the UNU.RAN library of standard distributions with truncated domain when the CDF is not available.
Get the area below the PDF of the distribution. If this area is not known,
unur_distr_cont_upd_pdfarea is called to compute it. If this is not successfulUNUR_INFINITYis returned andunur_errnois set toUNUR_ERR_DISTR_GET.
These are special cases of a continuous univariate distributions and thus they have most of these parameters (with the exception that functions cannot be changed). Additionally,
Create an object for order statistics of sample size n and rank k. distribution must be a pointer to a univariate continuous distribution. The resulting generator object is of the same type as of a unur_distr_cont_new call. (However, it cannot be used to make an order statistics out of an order statistics.)
To have a PDF for the order statistics, the given distribution object must contain a CDF and a PDF. Moreover, it is assumed that the given PDF is the derivative of the given CDF. Otherwise the area below the PDF of the order statistics is not computed correctly.
Important: There is no warning when the computed area below the PDF of the order statistics is wrong.
Get pointer to distribution object for underlying distribution.
Change sample size n and rank k of order statistics. In case of invalid data, no parameters are changed. The area below the PDF can be set to that of the underlying distribution by a unur_distr_corder_upd_pdfarea call.
Get sample size n and rank k of order statistics. In case of error an error code is returned.
Additionally most of the set and get calls for continuous univariate distributions work. The most important exceptions are that the PDF and CDF cannot be changed and unur_distr_cont_upd_mode uses in any way a (slow) numerical method that might fail.
— : UNUR_FUNCT_CONT* unur_distr_corder_get_dpdf (UNUR_DISTR* distribution)
— : UNUR_FUNCT_CONT* unur_distr_corder_get_cdf (UNUR_DISTR* distribution)
Get the respective pointer to the PDF, the derivative of the PDF and the CDF of the distribution, respectively. The pointer is of type
double funct(double x, UNUR_DISTR *distr). If the corresponding function is not available for the distribution, theNULLpointer is returned. See also unur_distr_cont_get_pdf. (Macro)
— : double unur_distr_corder_eval_dpdf (double x, UNUR_DISTR* distribution)
— : double unur_distr_corder_eval_cdf (double x, UNUR_DISTR* distribution)
Evaluate the PDF, derivative of the PDF. and the CDF, respectively, at x. Notice that distribution must not be the
NULLpointer. If the corresponding function is not available for the distribution,UNUR_INFINITYis returned andunur_errnois set toUNUR_ERR_DISTR_DATA. See also unur_distr_cont_eval_pdf. (Macro)IMPORTANT: In the case of a truncated standard distribution these calls always return the respective values of the untruncated distribution!
Set array of parameters for underlying distribution. See unur_distr_cont_set_pdfparams for details. (Macro)
Get number of parameters of the PDF of the underlying distribution and set pointer params to array of parameters. See unur_distr_cont_get_pdfparams for details. (Macro)
Set the left and right borders of the domain of the distribution. See unur_distr_cont_set_domain for details. (Macro)
Get the left and right borders of the domain of the distribution. See unur_distr_cont_get_domain for details. (Macro)
Get the left and right borders of the (truncated) domain of the distribution. See unur_distr_cont_get_truncated for details. (Macro)
The following paramters must be set whenever one of the essential parameters has been set or changed (and the parameter is required for the chosen method).
Set mode of distribution. See also unur_distr_corder_set_mode. (Macro)
Recompute the mode of the distribution numerically. Notice that this routine is slow and might not work properly in every case. See also unur_distr_cont_upd_mode for further details. (Macro)
Get mode of distribution. See unur_distr_cont_get_mode for details. (Macro)
Set the area below the PDF. See unur_distr_cont_set_pdfarea for details. (Macro)
Recompute the area below the PDF of the distribution. It only works for order statistics for distribution objects from the UNU.RAN library of standard distributions when the corresponding function is available. unur_distr_cont_upd_pdfarea assumes that the PDF of the underlying distribution is normalized, i.e. it is the derivative of its CDF. Otherwise the computed area is wrong and there is no warning about this failure. See unur_distr_cont_upd_pdfarea for further details. (Macro)
Get the area below the PDF of the distribution. See unur_distr_cont_get_pdfarea for details. (Macro)
Empirical univariate distributions are derived from observed data. There are two ways to create such a generator object:
Create a new (empty) object for empirical univariate continuous distribution.
Set observed sample for empirical distribution.
Read data from file filename. It reads the first number from each line. Numbers are parsed by means of the C standard routine strtod. Lines that do not start with
+,-,., or a digit are ignored. (Beware of lines starting with a blank!)In case of an error (file cannot be opened, invalid string for double in line) no data are copied into the distribution object and an error code is returned.
Get number of samples and set pointer sample to array of observations. If no sample has been given, an error code is returned and
sampleis set toNULL.Important: Do not change the entries in sample!
Set a histogram with bins of equal width. prob is an array of length n_prob that contains the probabilities for the bins (in ascending order). xmin and xmax give the lower and upper bound of the histogram, respectively. The bins are assumed to have equal width.
Remark: This is shortcut for calling unur_distr_cemp_set_hist_prob and unur_distr_cemp_set_hist_domain. Notice: All sampling methods either use raw data or histogram. It is possible to set both types of data; however, it is not checked whether the given histogran corresponds to possibly given raw data.
Set probabilities of a histogram with n_prob bins. Hence prob must be an array of length n_prob that contains the probabilities for the bins in ascending order. It is important also to set the location of the bins either with a unur_distr_cemp_set_hist_domain for bins of equal width or unur_distr_cemp_set_hist_bins when the bins have different width.
Notice: All sampling methods either use raw data or histogram. It is possible to set both types of data; however, it is not checked whether the given histogram corresponds to possibly given raw data.
Set a domain of a histogram with bins of equal width. xmin and xmax give the lower and upper bound of the histogram, respectively.
Set location of bins of a histogram with n_bins bins. Hence bins must be an array of length n_bins. The domain of the distribution is automatically set by this call and overrides any calls to unur_distr_cemp_set_hist_domain. Important: The probabilities of the bins of the distribution must be already be set by a unur_distr_cemp_set_hist_prob (or a unur_distr_cemp_set_hist call) and the value of n_bins must equal n_prob
+1from the corresponding value of the respective call.
The following calls handle multivariate distributions. However, the requirements of particular generation methods is not as unique as for univariate distributions. Moreover, random vector generation methods are still under development. The below functions are a first attempt to handle this situation.
Notice that some of the parameters – when given carelessly – might contradict to others. For example: Some methods require the marginal distribution and some methods need a standardized form of the marginal distributions, where the actual mean and variance is stored in the mean vector and the covariance matrix, respectively.
We also have to mention that some methods might abuse some of the parameters. Please read the discription of the chosen sampling method carfully.
The following kind of calls exists:
Some methods use the logarithm of the density if available.
Create a new (empty) object for multivariate continuous distribution. dim is the number of components of the random vector (i.e. its dimension). It is also possible to use dimension 1. Notice, however, that this is treated as a distribution of random vectors with only one component and not as a distribution of real numbers. For the latter unur_distr_cont_new should be used to create an object for a univariate distribution.
Set respective pointer to the PDF of the distribution. This function must be of type
double funct(const double *x, UNUR_DISTR *distr), where x must be a pointer to a double array of appropriate size (i.e. of the same size as given to the unur_distr_cvec_new call).It is not necessary that the given PDF is normalized, i.e. the integral need not be 1. Nevertheless the volume below the PDF can be provided by a unur_distr_cvec_set_pdfvol call.
It is not possible to change the PDF. Once the PDF is set it cannot be overwritten. This also holds when the logPDF is given. A new distribution object has to be used instead.
Set pointer to the gradient of the PDF. The type of this function must be
int funct(double *result, const double *x, UNUR_DISTR *distr), where result and x must be pointers to double arrays of appropriate size (i.e. of the same size as given to the unur_distr_cvec_new call). The gradient of the PDF is stored in the array result. The function should return an error code in case of an error and must returnUNUR_SUCCESSotherwise.The given function must be the gradient of the function given by a unur_distr_cvec_set_pdf call.
It is not possible to change the gradient of the PDF. Once the dPDF is set it cannot be overwritten. This also holds when the gradient of the logPDF is given. A new distribution object has to be used instead.
Set pointer to partial derivatives of the PDF. The type of this function must be
double funct(const double *x, int coord, UNUR_DISTR *distr), where x must be a pointer to a double array of appropriate size (i.e. of the same size as given to the unur_distr_cvec_new call). coord is the coordinate for which the partial dervative should be computed.Notice that coord must be an integer from {0,...,dim-1}.
It is not possible to change the partial derivative of the PDF. Once the pdPDF is set it cannot be overwritten. This also holds when the partial derivative of the logPDF is given. A new distribution object has to be used instead.
Get the pointer to the PDF of the distribution. The pointer is of type
double funct(const double *x, UNUR_DISTR *distr). If the corresponding function is not available for the distribution, theNULLpointer is returned.
Get the pointer to the gradient of the PDF of the distribution. The pointer is of type
int double funct(double *result, const double *x, UNUR_DISTR *distr). If the corresponding function is not available for the distribution, theNULLpointer is returned.
Evaluate the PDF of the distribution at x. x must be a pointer to a double array of appropriate size (i.e. of the same size as given to the unur_distr_cvec_new call) that contains the vector for which the function has to be evaluated.
Notice that distribution must not be the
NULLpointer. If the corresponding function is not available for the distribution,UNUR_INFINITYis returned andunur_errnois set toUNUR_ERR_DISTR_DATA.
Evaluate the gradient of the PDF of the distribution at x. The result is stored in the double array result. Both result and x must be pointer to double arrays of appropriate size (i.e. of the same size as given to the unur_distr_cvec_new call).
Notice that distribution must not be the
NULLpointer. If the corresponding function is not available for the distribution, an error code is returned andunur_errnois set toUNUR_ERR_DISTR_DATA(result is left unmodified).
Evaluate the partial derivative of the PDF of the distribution at x for the coordinate coord. x must be a pointer to a double array of appropriate size (i.e. of the same size as given to the unur_distr_cvec_new call) that contains the vector for which the function has to be evaluated.
Notice that coord must be an integer from {0,...,dim-1}.
Notice that distribution must not be the
NULLpointer. If the corresponding function is not available for the distribution,UNUR_INFINITYis returned andunur_errnois set toUNUR_ERR_DISTR_DATA.
— : int unur_distr_cvec_set_dlogpdf (UNUR_DISTR* distribution, UNUR_VFUNCT_CVEC* dlogpdf)
— : int unur_distr_cvec_set_pdlogpdf (UNUR_DISTR* distribution, UNUR_FUNCTD_CVEC* pdlogpdf)
— : UNUR_FUNCT_CVEC* unur_distr_cvec_get_logpdf (const UNUR_DISTR* distribution)
— : UNUR_VFUNCT_CVEC* unur_distr_cvec_get_dlogpdf (const UNUR_DISTR* distribution)
— : double unur_distr_cvec_eval_logpdf (const double* x, UNUR_DISTR* distribution)
— : int unur_distr_cvec_eval_dlogpdf (double* result, const double* x, UNUR_DISTR* distribution)
— : double unur_distr_cvec_eval_pdlogpdf (const double* x, int coord, UNUR_DISTR* distribution)
Analogous calls for the logarithm of the density function.
Set mean vector for multivariate distribution. mean must be a pointer to an array of size
dim, wheredimis the dimension returned by unur_distr_get_dim. ANULLpointer for mean is interpreted as the zero vector (0,...,0).Important: If the parameters of a distribution from the UNU.RAN library of standard distributions (see Standard distributions) are changed, then neither its mode nor the normalization constant are updated. Please use the respective calls unur_distr_cvec_upd_mode and unur_distr_cvec_upd_pdfvol.
Get the mean vector of the distribution. The function returns a pointer to an array of size
dim. If the mean vector is not marked as known theNULLpointer is returned andunur_errnois set toUNUR_ERR_DISTR_GET.Important: Do not modify the array that holds the mean vector!
Set covariance matrix for multivariate distribution. covar must be a pointer to an array of size
dimxdim, wheredimis the dimension returned by unur_distr_get_dim. The rows of the matrix have to be stored consecutively in this array.covar must be a variance-covariance matrix of the distribution, i.e. it must be symmetric and positive definit and its diagonal entries (i.e. the variance of the components of the random vector) must be strictly positive. The Cholesky factor is computed (and stored) to verify the positive definiteness condition. Notice that the inverse of the given covariance matrix is automatically computed when it is requested by some routine. Notice that the computation of this inverse matrix is unstable in case of high correlations and/or high dimensions. Thus it might fail and methods that require this inverse cannot be used. As an alternative the inverse of the covariance matrix can be directly set by a unur_distr_cvec_set_covar_inv call.
A
NULLpointer for covar is interpreted as the identity matrix.Important: This entry is abused in some methods which do not require the covariance matrix. It is then used to perform some transformation to obtain better performance.
Important: In case of an error (e.g. because covar is not a valid covariance matrix) an error code is returned. Moreover, the covariance matrix is not set and is marked as unknown. A previously set covariance matrix is then no longer available.
Important: If the parameters of a distribution from the UNU.RAN library of standard distributions (see Standard distributions) are changed, then neither its mode nor the normalization constant are updated. Please use the respective calls unur_distr_cvec_upd_mode and unur_distr_cvec_upd_pdfvol. Remark: UNU.RAN does not check whether the an eventually set covariance matrix and a rank-correlation matrix do not contradict each other.
Set inverse of the covariance matrix for multivariate distribution. covar_inv must be a pointer to an array of size
dimxdim, wheredimis the dimension returned by unur_distr_get_dim. The rows of the matrix have to be stored consecutively in this array.covar_inv must be symmetric and positive definit. Only the symmetry of the matrix is checked.
A
NULLpointer for covar_inv is interpreted as the identity matrix.Important: In case of an error (because covar_inv is not symetric) an error code is returned. Moreover, the inverse of the covariance matrix is not set and is marked as unknown. A previously set inverse matrix is then no longer available.
Remark: UNU.RAN does not check whether the given matrix is positive definit.
Remark: UNU.RAN does not check whether the matrix covar_inv is the inverse of the eventually set covariance matrix.
— : const double* unur_distr_cvec_get_cholesky (const UNUR_DISTR* distribution)
— : const double* unur_distr_cvec_get_covar_inv (UNUR_DISTR* distribution)
Get covariance matrix of distribution, its Cholesky factor, and its inverse, respectively. The function returns a pointer to an array of size
dimxdim. The rows of the matrix are stored consecutively in this array. If the requested matrix is not marked as known theNULLpointer is returned andunur_errnois set toUNUR_ERR_DISTR_GET.Important: Do not modify the array that holds the covariance matrix!
Remark: The inverse of the covariance matrix is computed if it is not already stored.
Set rank-correlation matrix (Spearman's correlation) for multivariate distribution. rankcorr must be a pointer to an array of size
dimxdim, wheredimis the dimension returned by unur_distr_get_dim. The rows of the matrix have to be stored consecutively in this array.rankcorr must be a rank-correlation matrix of the distribution, i.e. it must be symmetric and positive definite and its diagonal entries must be equal to
1.The Cholesky factor is computed (and stored) to verify the positive definiteness condition.
A
NULLpointer for rankcorr is interpreted as the identity matrix.Important: In case of an error (e.g. because rankcorr is not a valid rank-correlation matrix) an error code is returned. Moreover, the rank-correlation matrix is not set and is marked as unknown. A previously set rank-correlation matrix is then no longer available.
Remark: UNU.RAN does not check whether the an eventually set covariance matrix and a rank-correlation matrix do not contradict each other.
— : const double* unur_distr_cvec_get_rk_cholesky (const UNUR_DISTR* distribution)
Get rank-correlation matrix and its cholesky factor, respectively, of distribution. The function returns a pointer to an array of size
dimxdim. The rows of the matrix are stored consecutively in this array. If the requested matrix is not marked as known theNULLpointer is returned andunur_errnois set toUNUR_ERR_DISTR_GET.Important: Do not modify the array that holds the rank-correlation matrix!
Sets marginal distributions of the given distribution to the same marginal distribution object. The marginal distribution must be an instance of a continuous univariate distribution object. Notice that the marginal distribution is copied into the distribution object.
Analogously to the above unur_distr_cvec_set_marginals call. However, now an array marginals of the pointers to each of the marginal distributions must be given. It must be an array of size
dim, wheredimis the dimension returned by unur_distr_get_dim. Notice: Local copies for each of the entries are stored in the distribution object. If some of these entries are identical (i.e. contain the same pointer), then for each of these a new copy is made.
Similar to the above unur_distr_cvec_set_marginal_array call. However, now the pointers to the particular marginal distributions can be given as parameter and does not require an array of pointers. Additionally the given distribution objects are immediately destroyed. Thus calls like unur_distr_normal can be used as arguments. (With unur_distr_cvec_set_marginal_array the result of such call has to be stored in a pointer since it has to be freed afterwarts to avoid memory leaks!)
The number of pointers to in the list of function arguments must be equal to the dimension of the distribution, i.e. the dimension returned by unur_distr_get_dim. If one of the given pointer to marginal distributions is the
NULLpointer then the marginal distributions of distribution are not set (or previous settings are not changed) and an error code is returned.Important: All distribution objects given in the argument list are destroyed!
Get pointer to the n-th marginal distribution object from the given multivariate distribution. If this does not exist,
NULLis returned. The marginal distributions are enumerated from1todim, wheredimis the dimension returned by unur_distr_get_dim.
Sets array of parameters for distribution. There is an upper limit for the number of parameters
n_params. It is given by the macroUNUR_DISTR_MAXPARAMSin unuran_config.h. (It is set to 5 by default but can be changed to any appropriate nonnegative number.) If n_params is negative or exceeds this limit no parameters are copied into the distribution object andunur_errnois set toUNUR_ERR_DISTR_NPARAMS.For standard distributions from the UNU.RAN library the parameters are checked. Moreover, the domain is updated automatically. If the given parameters are invalid for the standard distribution, then no parameters are set and an error code is returned. Notice that the given parameter list for such a distribution is handled in the same way as in the corresponding new calls, i.e. optional parameters for the PDF that are not present in the given list are (re-)set to their default values.
Important: If the parameters of a distribution from the UNU.RAN library of standard distributions (see Standard distributions) are changed, then neither its mode nor the normalization constant are updated. Please use the respective calls unur_distr_cvec_upd_mode and unur_distr_cvec_upd_pdfvol.
Get number of parameters of the PDF and set pointer params to array of parameters. If no parameters are stored in the object, an error code is returned and
paramsis set toNULL.Important: Do not change the entries in params!
This function provides an interface for additional vector parameters for a multivariate distribution besides mean vector and covariance matrix which have their own calls.
It sets the parameter with number par. par indicates directly which of the parameters is set and must be a number between
0andUNUR_DISTR_MAXPARAMS-1 (the upper limit of possible parameters defined in unuran_config.h; it is set to 5 but can be changed to any appropriate nonnegative number.)The entries of a this parameter are given by the array param_vec of size n_params. Notice that using this interface an An (n x m)-matrix has to be stored in an array of length n_params = n times m; where the rows of the matrix are stored consecutively in this array.
Due to great variety of possible parameters for a multivariate distribution there is no simpler interface.
If param_vec is
NULLthen the corresponding entry is cleared.Important: If the parameters of a distribution from the UNU.RAN library of standard distributions (see Standard distributions) are changed, then neither its mode nor the normalization constant are updated. Please use the respective calls unur_distr_cvec_upd_mode and unur_distr_cvec_upd_pdfvol. If an error occurs no parameters are copied into the parameter object
unur_errnois set toUNUR_ERR_DISTR_DATA.
Get parameter of the PDF with number par. The pointer to the parameter array is stored in param_vecs, its size is returned by the function. If the requested parameter is not set, then an error code is returned and
paramsis set toNULL.Important: Do not change the entries in param_vecs!
Set rectangular domain for distribution with lowerleft and upperright vertices. Both must be pointer to an array of the size returned by unur_distr_get_dim. A
NULLpointer is interpreted as the zero vector (0,...,0). For setting a coordinate of the boundary to +/- infinityuse+/- UNUR_INFINITY. The lowerleft vertex must be strictly smaller than upperright in each component. Otherwise no domain is set andunur_errnois set toUNUR_ERR_DISTR_SET.By default the domain of a distribution is unbounded. Thus one can use this call to truncate an existing distribution.
Important: Changing the domain of distribution marks derived parameters like the mode or the center as unknown and must be set after changing the domain. This is important for the already set (or default) value for the center does not fall into the given domain. Notice that calls of the PDF and derived functions return
0.when the parameter is not contained in the domain.
Check whether x falls into the domain of distribution.
The following paramters must be set whenever one of the essential parameters has been set or changed (and the parameter is required for the chosen method).
Set mode of the distribution. mode must be a pointer to an array of the size returned by unur_distr_get_dim. A
NULLpointer for mode is interpreted as the zero vector (0,...,0).
Recompute the mode of the distribution. This call works properly for distribution objects from the UNU.RAN library of standard distributions when the corresponding function is available. If it failes
unur_errnois set toUNUR_ERR_DISTR_DATA.
Get mode of the distribution. The function returns a pointer to an array of the size returned by unur_distr_get_dim. If the mode is not marked as known the
NULLpointer is returned andunur_errnois set toUNUR_ERR_DISTR_GET. (There is no difference between the case where no routine for computing the mode is available and the case where no mode exists for the distribution at all.)Important: Do not modify the array that holds the mode!
Set center of the distribution. center must be a pointer to an array of the size returned by unur_distr_get_dim. A
NULLpointer for center is interpreted as the zero vector (0,...,0).The center is used by some methods to shift the distribution in order to decrease numerical round-off error. If not given explicitly a default is used. Moreover, it is used as starting point for several numerical search algorithm (e.g. for the mode). Then center must be a pointer where the call to the PDF returns a non-zero value. In particular center must contained in the domain of the distribution.
Default: The mode, if given by a unur_distr_cvec_set_mode call; else the mean, if given by a unur_distr_cvec_set_mean call; otherwise the null vector (0,...,0).
Get center of the distribution. The function returns a pointer to an array of the size returned by unur_distr_get_dim. It always returns some point as there always exists a default for the center, see unur_distr_cvec_set_center. Important: Do not modify the array that holds the center!
Set the volume below the PDF. If vol is non-positive, no volume is set and
unur_errnois set toUNUR_ERR_DISTR_SET.
Recompute the volume below the PDF of the distribution. It only works for distribution objects from the UNU.RAN library of standard distributions when the corresponding function is available. Otherwise
unur_errnois set toUNUR_ERR_DISTR_DATA.This call also sets the normalization constant such that the given PDF is the derivative of a given CDF, i.e. the volume is 1.
Get the volume below the PDF of the distribution. If this volume is not known,
unur_distr_cont_upd_pdfarea is called to compute it. If this is not successfulUNUR_INFINITYis returned andunur_errnois set toUNUR_ERR_DISTR_GET.
Full conditional distribution for a given continuous multivariate distributiion. The condition is a position vector and either a variable that is variated or a vector that indicates the direction on which the random vector can variate.
There is a subtle difference between using direction vector and using the k-th variable. When a direction vector is given the PDF of the conditional distribution is defined by f(t) = PDF(pos + t * dir). When a variable is selected the full conditional distribution with all other variables fixed is used.
This is a special case of a continuous univariate distribution and thus they have most of these parameters (with the exception that functions cannot be changed). Additionally,
This distibution type is primarily used for evaluation the conditional distribution and its derivative (as required for, e.g., the Gibbs sampler). The density is not normalized (i.e. does not integrate to one). Mode and area are not available and it does not make sense to use any call to set or change parameters except the ones given below.
Create an object for full conditional distribution for the given distribution. The condition is given by a position vector pos and either the k-th variable that is variated or the vector dir that contains the direction on which the random vector can variate.
distribution must be a pointer to a multivariate continuous distribution. pos must be a pointer to an array of size
dim, wheredimis the dimension of the underlying distribution object. dir must be a pointer to an array if sizedimorNULL. k must be in the range0, ..., dim-1. If the k-th variable is used, dir must be set toNULL.Notice: There is a subtle difference between using direction vector dir and using the k-th variable. When dir is given, the current position pos is mapped into 0 of the conditional distribution and the derivative is taken from the function PDF(pos+t*dir) w.r.t. t. On the other hand, when the coordinate k is used (i.e., when dir is set to
NULL), the full conditional distribution of the distribution is considered (as used for the Gibbs sampler). In particular, the current point is just projected into the one-dimensional subspace without mapping it into the point 0.Notice: If a coordinate k is used, then the k-th partial derivative is used if it as available. Otherwise the gradient is computed and the k-th component is returned.
The resulting generator object is of the same type as of a unur_distr_cont_new call.
Set/change condition for conditional distribution. Change values of fixed variables to pos and use direction dir or k-th variable of conditional distribution.
pos must be a pointer to an array of size
dim, wheredimis the dimension of the underlying distribution object. dir must be a pointer to an array if sizedimorNULL. k must be in the range0, ..., dim-1. If the k-th variable is used, dir must be set toNULL.Notice: There is a subtle difference between using direction vector dir and using the k-th variable. When dir is given, the current position pos is mapped into 0 of the conditional distribution and the derivative is taken from the function PDF(pos+t*dir) w.r.t. t. On the other hand, when the coordinate k is used (i.e., when dir is set to
NULL), the full conditional distribution of the distribution is considered (as used for the Gibbs sampler). In particular, the current point is just projected into the one-dimensional subspace without mapping it into the point 0.
Get condition for conditional distribution. The values for the fixed variables are stored in pos, which must be a pointer to an array of size
dim. The condition is stored in dir and k, respectively.Important: Do not change the entries in pos and dir!
Get pointer to distribution object for underlying distribution.
Empirical multivariate distributions are just lists of vectors (with the same dimension). Thus there are only calls to insert these data. How these data are used to sample from the empirical distribution depends from the chosen generation method.
Create a new (empty) object for an empirical multivariate continuous distribution. dim is the number of components of the random vector (i.e. its dimension). It must be at least 2; otherwise unur_distr_cemp_new should be used to create an object for an empirical univariate distribution.
Set observed sample for empirical distribution. sample is an array of doubles of size
dimx n_sample, wheredimis the dimension of the distribution returned by unur_distr_get_dim. The data points must be stored consecutively in sample, i.e., data points (x1, y1), (x2, y2), ... are given as an array {x1, y1, x2, y2, ...}.
Read data from file filename. It reads the first
dimnumbers from each line, wheredimis the dimension of the distribution returned by unur_distr_get_dim. Numbers are parsed by means of the C standard routine strtod. Lines that do not start with+,-,., or a digit are ignored. (Beware of lines starting with a blank!)In case of an error (file cannot be opened, too few entries in a line, invalid string for double in line) no data are copied into the distribution object and an error code is returned.
Get number of samples and set pointer sample to array of observations. If no sample has been given, an error code is returned and sample is set to
NULL. If successful sample points to an array of lengthdimxn_sample, wheredimis the dimension of the distribution returned by unur_distr_get_dim andn_samplethe return value of the function.Important: Do not modify the array sample.
Distributions for random matrices. Notice that UNU.RAN uses
arrays of doubles to handle matrices. The rows of
the matrix are stored consecutively.
Create a new (empty) object for a matrix distribution. n_rows and n_cols are the respective numbers of rows and columns of the random matrix (i.e. its dimensions). It is also possible to have only one number or rows and/or columns. Notice, however, that this is treated as a distribution of random matrices with only one row or column or component and not as a distribution of vectors or real numbers. For the latter unur_distr_cont_new or unur_distr_cvec_new should be used to create an object for a univariate distribution and a multivariate (vector) distribution, respectively.
Get number of rows and columns of random matrix (its dimension). It returns the total number of components. If successfull
UNUR_SUCCESSis returned.
The calls in this section can be applied to discrete univariate distributions.
doubles.
It can be automatically computed if the pmf is
given but pv is not.
Create a new (empty) object for a univariate discrete distribution.
There are two interfaces for discrete univariate distributions: Either provide a (finite) probability vector (PV). Or provide a probability mass function (PMF). For the latter case there are also a couple of derived parameters that are not required when a PV is given.
It is not possible to set both a PMF and a PV directly. However, the PV can be computed from the PMF (or the CDF if no PMF is available) by means of a unur_distr_discr_make_pv call. If both the PV and the PMF are given in the distribution object it depends on the generation method which of these is used.
Set finite probability vector (PV) for the distribution. It is not necessary that the entries in the given PV sum to 1. n_pv must be positive. However, there is no testing whether all entries in pv are non-negative.
If no domain has been set, then the left boundary is set to
0, by default. If n_pv is too large, e.g. because left boundary + n_pv exceeds the range of integers, then the call fails.Notice that it is not possible to set both a PV and a PMF or CDF. If the PMF or CDF is set first one cannot set the PV. If the PMF or CDF is set first after a PV is set, the latter is removed (and recomputed using unur_distr_discr_make_pv when required).
Compute a PV when a PMF or CDF is given. However, when the domain is not given or is too large and the sum over the PMF is given then the (right) tail of the distribution is chopped off such that the probability for the tail region is less than 1.e-8. If the sum over the PMF is not given a PV of maximal length is computed.
The maximal size of the created PV is bounded by the macro
UNUR_MAX_AUTO_PVthat is defined in unuran_config.h.If successful, the length of the generated PV is returned. If the sum over the PMF on the chopped tail is not neglible small (i.e. greater than 1.e-8 or unknown) than the negative of the length of the PV is returned and
unur_errnois set toUNUR_ERR_DISTR_SET.Notice that the left boundary of the PV is set to
0by default when a discrete distribution object is created from scratch.If computing a PV fails for some reasons, an error code is returned and
unur_errnois set toUNUR_ERR_DISTR_SET.
Get length of PV of the distribution and set pointer pv to array of probabilities. If no PV is given, an error code is returned and pv is set to
NULL.
(It does not call unur_distr_discr_make_pv !)
— : int unur_distr_discr_set_cdf (UNUR_DISTR* distribution, UNUR_FUNCT_DISCR* cdf)
Set respective pointer to the PMF and the CDF of the distribution. These functions must be of type
double funct(int k, const UNUR_DISTR *distr).It is important to note that all these functions must return a result for all integers k. E.g., if the domain of a given PMF is the interval {1,2,3,...,100}, than the given function must return
0.0for all points outside this interval.The default domain for the PMF or CDF is [
0,INT_MAX]. The domain can be changed using a unur_distr_discr_set_domain call.It is not possible to change such a function. Once the PMF or CDF is set it cannot be overwritten. A new distribution object has to be used instead.
Notice that it is not possible to set both a PV and a PMF or CDF. If the PMF or CDF is set first one cannot set the PV. If the PMF or CDF is set first after a PV is set, the latter is removed (and recomputed using unur_distr_discr_make_pv when required).
— : double unur_distr_discr_eval_pmf (int k, const UNUR_DISTR* distribution)
— : double unur_distr_discr_eval_cdf (int k, const UNUR_DISTR* distribution)
Evaluate the PV, PMF, and the CDF, respectively, at k. Notice that distribution must not be the
NULLpointer. If no PV is set for the distribution, then unur_distr_discr_eval_pv behaves like unur_distr_discr_eval_pmf. If the corresponding function is not available for the distribution,UNUR_INFINITYis returned andunur_errnois set toUNUR_ERR_DISTR_DATA.IMPORTANT: In the case of a truncated standard distribution these calls always return the respective values of the untruncated distribution!
This function provides an alternative way to set a PMF of the distribution. pmfstr is a character string that contains the formula for the PMF, see Function String, for details. See also the remarks for the unur_distr_discr_set_pmf call.
It is not possible to call this funtion twice or to call this function after a unur_distr_discr_set_pmf call.
This function provides an alternative way to set a CDF; analogously to the unur_distr_discr_set_pmfstr call.
— : char* unur_distr_discr_get_cdfstr (const UNUR_DISTR* distribution)
Get pointer to respective string for PMF and CDF of distribution that is given via the string interface. This call allocates memory to produce this string. It should be freed when it is not used any more.
Set array of parameters for distribution. There is an upper limit for the number of parameters n_params. It is given by the macro
UNUR_DISTR_MAXPARAMSin unuran_config.h. (It is set to 5 but can be changed to any appropriate nonnegative number.) If n_params is negative or exceeds this limit no parameters are copied into the distribution object andunur_errnois set toUNUR_ERR_DISTR_NPARAMS.For standard distributions from the UNU.RAN library the parameters are checked. Moreover, the domain is updated automatically unless it has been changed before by a unur_distr_discr_set_domain call. If the given parameters are invalid for the standard distribution, then no parameters are set and an error code is returned. Notice that the given parameter list for such a distribution is handled in the same way as in the corresponding new calls, i.e. optional parameters for the PDF that are not present in the given list are (re-)set to their default values.
Important: Integer parameter must be given as
doubles.
Get number of parameters of the PMF and set pointer params to array of parameters. If no parameters are stored in the object, an error code is returned and
paramsis set toNULL.
Set the left and right borders of the domain of the distribution. This can also be used to truncate an existing distribution. For setting the boundary to +/- infinityuse
INT_MINandINT_MAX, respectively. If right is not strictly greater than left no domain is set andunur_errnois set toUNUR_ERR_DISTR_SET. It is allowed to use this call to increase the domain. If the PV of the discrete distribution is used, than the right boudary is ignored (and internally set to left + size of PV - 1). Notice thatINT_MINandINT_MAXare interpreted as (minus/plus) infinity.Default: [
0,INT_MAX].
Get the left and right borders of the domain of the distribution. If the domain is not set explicitly the interval [
INT_MIN,INT_MAX] is assumed and returned. When a PV is given then the domain is set automatically to [0,size of PV - 1].
The following paramters must be set whenever one of the essential parameters has been set or changed (and the parameter is required for the chosen method).
Recompute the mode of the distribution. This call works properly for distribution objects from the UNU.RAN library of standard distributions when the corresponding function is available. Otherwise a (slow) numerical mode finder is used. It only works properly for unimodal probability mass functions. If it failes
unur_errnois set toUNUR_ERR_DISTR_DATA.
Get mode of distribution. If the mode is not marked as known, unur_distr_discr_upd_mode is called to compute the mode. If this is not successful
INT_MAXis returned andunur_errnois set toUNUR_ERR_DISTR_GET. (There is no difference between the case where no routine for computing the mode is available and the case where no mode exists for the distribution at all.)
Set the sum over the PMF. If
sumis non-positive, no sum is set andunur_errnois set toUNUR_ERR_DISTR_SET.For a distribution object created by the UNU.RAN library of standard distributions you always should use the unur_distr_discr_upd_pmfsum. Otherwise there might be ambiguous side-effects.
Recompute the sum over the PMF of the distribution. In most cases the normalization constant is recomputed and thus the sum is 1. This call works for distribution objects from the UNU.RAN library of standard distributions when the corresponding function is available. When a PV, a PMF with finite domain, or a CDF is given, a simple generic function which uses a naive summation loop is used. If this computation is not possible, an error code is returned and
unur_errnois set toUNUR_ERR_DISTR_DATA.The call does not work for distributions from the UNU.RAN library of standard distributions with truncated domain when the CDF is not available.
Get the sum over the PMF of the distribution. If this sum is not known, unur_distr_discr_upd_pmfsum is called to compute it. If this is not successful
UNUR_INFINITYis returned andunur_errnois set toUNUR_ERR_DISTR_GET.
Sampling from a particular distribution with UNU.RAN requires the following steps:
Important: Initialization of the generator object might fail.
unur_init
returns a NULL pointer then, which must not be
used for sampling.
Important: Currently not all methods allow reinitialization, see the description of the particular method (keyword Reinit).
Important: Reinitialization of the generator object might fail. Thus one must check the return code of the unur_reinit call.
Important: When reinitialization fails then sampling routines
always return INFINITY (for continuous distributions) or 0 (for
discrete distributions), respectively.
However, it is still possible to change the underlying distribution
and try to reinitialize again.
Routines for all generator objects.
Initialize a generator object. All necessary information must be stored in the parameter object.
Important: If an error has occurred a
NULLpointer is return. This must not be used for the sampling routines (this causes a segmentation fault).Always check whether the call was successful or not!
Important: This call destroys the parameter object automatically. Thus it is not necessary/allowed to free it.
Update an existing generator object after the underlying distribution has been modified (using unur_get_distr together with corresponding set calls. It must be executed before sampling using this generator object is continued as otherwise it produces an invalid sample or might even cause a segmentation fault.
Important: Currently not all methods allow reinitialization, see the description of the particular method (keyword Reinit).
Important: Reinitialization of the generator object might fail. Thus one must check the return code:
UNUR_SUCCESS (0x0u)- success (no error)
UNUR_ERR_NO_REINIT- reinit routine not implemented.
- other values
- some error has occured while trying to reinitialize the generator object.
Important: When reinitialization fails then sampling routines always return INFINITY (for continuous distributions) or
0(for discrete distributions), respectively. However, it is still possible to change the underlying distribution and try to reinitialize again.Important: When one tries to run unur_reinit but reinitialization is not implemented, then the generator object cannot be used any more and must be destroyed and a new one has to be built from scratch.
— : double unur_sample_cont (UNUR_GEN* generator)
— : int unur_sample_vec (UNUR_GEN* generator, double* vector)
— : int unur_sample_matr (UNUR_GEN* generator, double* matrix)
Sample from generator object. The three routines depend on the type of the generator object (discrete or continuous univariate distribution, multivariate distribution, or random matrix).
Notice: UNU.RAN uses arrays of
doubles to handle matrices. There the rows of the matrix are stored consecutively.Notice: The routines unur_sample_vec and unur_sample_matr return
UNUR_SUCCESSif generation was successful and some error code otherwise.Important: These routines do not check whether generator is an invalid
NULLpointer.
Compute the U quantile of a continuous distribution using a generator object that implements an (approximate) inversion methods.
The following methods are currently available:
- HINV, see HINV.
- NINV, see NINV.
- PINV, see PINV.
- CSTD, see CSTD.
This requires that generator implements an inversion method.- DGT, see DGT.
The return value is (of course) type casted todouble.Important: This routine does not check whether generator is an invalid
NULLpointer.In case of an error UNUR_INFINITY or INT_MAX (depending on the type of generator) is returned.
Get a string with informations about the given generator. These informations allow some fine tuning of the generation method. If help is
TRUE, some hints on setting parameters are given.This function is intented for using in interactive environments (like R).
If an error occurs, then
NULLis returned.
Get the number of dimension of a (multivariate) distribution. For a univariate distribution
1is return.
Get pointer to distribution object from generator object. This function can be used to change the parameters of the distribution and reinitialize the generator object. Notice that currently not all generating methods have a reinitialize routine. This function should be used with extreme care. Changing the distribution is changed and using the generator object without reinitializing might cause wrong samples or segmentation faults. Moreover, if the corresponding generator object is freed, the pointer must not be used.
Important: The returned distribution object must not be freed. If the distribution object is changed then one must run unur_reinit !
Set flag whether the generator object should make a private copy of the given distribution object or just stores the pointer to this distribution object. Values for use_privatecopy:
TRUE- make a private copy (default)
FALSE- do not make a private copy and store pointer to given (external) distribution object.
By default, generator objects keep their own private copy of the given distribution object. Thus the generator object can be handled independently from other UNU.RAN objects (with uniform random number generators as the only exception). When the generator object is initialized the given distribution object is cloned and stored.
However, in some rare situations it can be useful when only the pointer to the given distribution object is stored without making a private copy. A possible example is when only one random variate has to be drawn from the distribution. This behavior can be achieved when use_localcopy is set to
FALSE.Warning! Using a pointer to the external distribution object instead of a private copy must be done with extreme care! When the distrubtion object is changed or freed then the generator object does not work any more, might case a segmentation fault, or (even worse) produces garbage. On the other hand, when the generator object is initialized or used to draw a random sampling the distribution object may be changed.
Notice: The prototypes of all
unur_<method>_newcalls use aconstqualifier for the distribution argument. However, if use_privatecopy is set toFALSEthis qualifier is discarded and the distribution might be changed.Important! If use_localcopy is set to
FALSEand the corresponding distribution object is changed then one must run unur_reinit on the generator object. (Notice that currently not all generation methods support reinitialization.)Default: use_privatecopy is
TRUE.
AUTO selects a an appropriate method for the given distribution object automatically. There are no parameters for this method, yet. But it is planned to give some parameter to describe the task for which the random variate generator is used for and thus make the choice of the generating method more appropriate. Notice that the required sampling routine for the generator object depends on the type of the given distribution object.
The chosen method also depends on the sample size for which the generator object will be used. If only a few random variates the order of magnitude of the sample size should be set via a unur_auto_set_logss call.
IMPORTANT: This is an experimental version and the method chosen may change in future releases of UNU.RAN.
For an example see As short as possible.
Create a generator object for the given distribution object.
Set the order of magnitude for the size of the sample that will be generated by the generator, i.e., the the common logarithm of the sample size.
Default is 10.
Notice: This feature will be used in future releases of UNU.RAN only.
|
Methods for continuous univariate distributions sample with unur_sample_cont
|
/* ------------------------------------------------------------- */
/* File: example_cont.c */
/* ------------------------------------------------------------- */
/* Include UNURAN header file. */
#include <unuran.h>
/* ------------------------------------------------------------- */
/* Example how to sample from a continuous univariate */
/* distribution. */
/* */
/* We build a distribution object from scratch and sample. */
/* ------------------------------------------------------------- */
/* Define the PDF and dPDF of our distribution. */
/* */
/* Our distribution has the PDF */
/* */
/* / 1 - x*x if |x| <= 1 */
/* f(x) = < */
/* \ 0 otherwise */
/* */
/* The PDF of our distribution: */
double mypdf( double x, const UNUR_DISTR *distr )
/* The second argument (`distr') can be used for parameters */
/* for the PDF. (We do not use parameters in our example.) */
{
if (fabs(x) >= 1.)
return 0.;
else
return (1.-x*x);
} /* end of mypdf() */
/* The derivative of the PDF of our distribution: */
double mydpdf( double x, const UNUR_DISTR *distr )
{
if (fabs(x) >= 1.)
return 0.;
else
return (-2.*x);
} /* end of mydpdf() */
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double x; /* will hold the random number */
/* Declare the three UNURAN objects. */
UNUR_DISTR *distr; /* distribution object */
UNUR_PAR *par; /* parameter object */
UNUR_GEN *gen; /* generator object */
/* Create a new distribution object from scratch. */
/* Get empty distribution object for a continuous distribution */
distr = unur_distr_cont_new();
/* Fill the distribution object -- the provided information */
/* must fulfill the requirements of the method choosen below. */
unur_distr_cont_set_pdf(distr, mypdf); /* PDF */
unur_distr_cont_set_dpdf(distr, mydpdf); /* its derivative */
unur_distr_cont_set_mode(distr, 0.); /* mode */
unur_distr_cont_set_domain(distr, -1., 1.); /* domain */
/* Choose a method: TDR. */
par = unur_tdr_new(distr);
/* Set some parameters of the method TDR. */
unur_tdr_set_variant_gw(par);
unur_tdr_set_max_sqhratio(par, 0.90);
unur_tdr_set_c(par, -0.5);
unur_tdr_set_max_intervals(par, 100);
unur_tdr_set_cpoints(par, 10, NULL);
/* Create the generator object. */
gen = unur_init(par);
/* Notice that this call has also destroyed the parameter */
/* object `par' as a side effect. */
/* It is important to check if the creation of the generator */
/* object was successful. Otherwise `gen' is the NULL pointer */
/* and would cause a segmentation fault if used for sampling. */
if (gen == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* It is possible to reuse the distribution object to create */
/* another generator object. If you do not need it any more, */
/* it should be destroyed to free memory. */
unur_distr_free(distr);
/* Now you can use the generator object `gen' to sample from */
/* the distribution. Eg.: */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* When you do not need the generator object any more, you */
/* can destroy it. */
unur_free(gen);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
/* ------------------------------------------------------------- */
/* File: example_cont_str.c */
/* ------------------------------------------------------------- */
/* String API. */
/* ------------------------------------------------------------- */
/* Include UNURAN header file. */
#include <unuran.h>
/* ------------------------------------------------------------- */
/* Example how to sample from a continuous univariate */
/* distribution. */
/* We use a generic distribution object and sample. */
/* */
/* The PDF of our distribution is given by */
/* */
/* / 1 - x*x if |x| <= 1 */
/* f(x) = < */
/* \ 0 otherwise */
/* */
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double x; /* will hold the random number */
/* Declare UNURAN generator object. */
UNUR_GEN *gen; /* generator object */
/* Create the generator object. */
/* Use a generic continuous distribution. */
/* Choose a method: TDR. */
gen = unur_str2gen(
"distr = cont; pdf=\"1-x*x\"; domain=(-1,1); mode=0. & \
method=tdr; variant_gw; max_sqhratio=0.90; c=-0.5; \
max_intervals=100; cpoints=10");
/* It is important to check if the creation of the generator */
/* object was successful. Otherwise `gen' is the NULL pointer */
/* and would cause a segmentation fault if used for sampling. */
if (gen == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* Now you can use the generator object `gen' to sample from */
/* the distribution. Eg.: */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* When you do not need the generator object any more, you */
/* can destroy it. */
unur_free(gen);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
AROU is a variant of the ratio-of-uniforms method that uses the fact that the transformed region is convex for many distributions. It works for all T-concave distributions with T(x) = -1/sqrt(x).
It is possible to use this method for correlation induction by setting an auxiliary uniform random number generator via the unur_set_urng_aux call. (Notice that this must be done after a possible unur_set_urng call.) When an auxiliary generator is used then the number of used uniform random numbers that is used up for one generated random variate is constant and equal to 1.
There exists a test mode that verifies whether the conditions for the method are satisfied or not while sampling. It can be switched on by calling unur_arou_set_verify and unur_arou_chg_verify respectively. Notice however that sampling is (much) slower then.
For densities with modes not close to 0 it is suggested to set either the mode or the center of the distribution by the unur_distr_cont_set_mode or unur_distr_cont_set_center call. The latter is the approximate location of the mode or the mean of the distribution. This location provides some information about the main part of the PDF and is used to avoid numerical problems.
If usedars is set to
TRUE, “derandomized adaptive rejection sampling” (DARS) is used in setup. Segments where the area between hat and squeeze is too large compared to the average area between hat and squeeze over all intervals are split. This procedure is repeated until the ratio between area below squeeze and area below hat exceeds the bound given by unur_arou_set_max_sqhratio call or the maximum number of segments is reached. Moreover, it also aborts when no more segments can be found for splitting.Segments are split such that the angle of the segments are halved (corresponds to arc-mean rule of method TDR (see TDR)).
Default is
TRUE.
Set factor for “derandomized adaptive rejection sampling”. This factor is used to determine the segments that are “too large”, that is, all segments where the area between squeeze and hat is larger than factor times the average area over all intervals between squeeze and hat. Notice that all segments are split when factor is set to
0., and that there is no splitting at all when factor is set toUNUR_INFINITY.Default is
0.99. There is no need to change this parameter.
Set upper bound for the ratio (area inside squeeze) / (area inside envelope). It must be a number between 0 and 1. When the ratio exceeds the given number no further construction points are inserted via adaptive rejection sampling. Use
0if no construction points should be added after the setup. Use1if adding new construction points should not be stopped until the maximum number of construction points is reached.Default is
0.99.
Get the current ratio (area inside squeeze) / (area inside envelope) for the generator. (In case of an error
UNUR_INFINITYis returned.)
Get the area below the hat for the generator. (In case of an error
UNUR_INFINITYis returned.)
Get the area below the squeeze for the generator. (In case of an error
UNUR_INFINITYis returned.)
Set maximum number of segements. No construction points are added after the setup when the number of segments succeeds max_segs.
Default is
100.
Set construction points for enveloping polygon. If stp is
NULL, then a heuristical rule of thumb is used to get n_stp construction points. This is the default behavior when this routine is not called. The (default) number of construction points is30, then.
Use the center as construction point. Default is
TRUE.
Set factor for relative size of the guide table for indexed search (see also method DGT DGT). It must be greater than or equal to
0. When set to0, then sequential search is used.Default is
2.
— : int unur_arou_chg_verify (UNUR_GEN* generator, int verify)
Turn verifying of algorithm while sampling on/off. If the condition squeeze(x) <= PDF(x) <= hat(x) is violated for some x then
unur_errnois set toUNUR_ERR_GEN_CONDITION. However notice that this might happen due to round-off errors for a few values of x (less than 1%).Default is
FALSE.
Sometimes it might happen that unur_init has been executed successfully. But when additional construction points are added by adaptive rejection sampling, the algorithm detects that the PDF is not T-concave.
With pedantic being
TRUE, the sampling routine is then exchanged by a routine that simply returnsUNUR_INFINITY. Otherwise the new point is not added to the list of construction points. At least the hat function remains T-concave.Setting pedantic to
FALSEallows sampling from a distribution which is “almost” T-concave and small errors are tolerated. However it might happen that the hat function cannot be improved significantly. When the hat function that has been constructed by the unur_init call is extremely large then it might happen that the generation times are extremely high (even hours are possible in extremely rare cases).Default is
FALSE.
ARS is an acceptance/rejection method that uses the concavity
of the log-density function to construct hat function and
squeezes automatically.
It is very similar to method TDR (see TDR) with variant GW,
parameter c = 0, and DARS switched off.
Moreover, method ARS requires the logPDF and its derivative
dlogPDF to run. On the other hand, it is designed to draw only a
(very) small samples and it is much more robust against
densities with very large or small areas below the PDF as
it occurs, for example, in conditional distributions of
(high dimensional) multivariate distributions.
Additionally, it can be re-initialized when the underlying
distribution has been modified.
Thus it is well suited for Gibbs sampling.
Notice, that method ARS is a restricted version of TDR. If the full functionally of Transformed Density Rejection is needed use method TDR.
Method ARS is designed for distributions with log-concave densities. To use this method you need a distribution object with the logarithm of the PDF and its derivative given.
The number of construction points as well as a set of such points can be provided using unur_ars_set_cpoints. Notice that addition construction points are added by means of adaptive rejection sampling until the maximal number of intervals given by unur_ars_set_max_intervals is reached.
A generated distribution object can be reinitialized using the unur_reinit call. When unur_reinit is called construction points for the new generator are necessary. There are two options: Either the same construction points as for the initial generator (given by a unur_ars_set_cpoints call) are used (this is the default), or percentiles of the old hat function can be used. This can be set or changed using unur_ars_set_reinit_percentiles and unur_ars_chg_reinit_percentiles. This feature is usefull when the underlying distribution object is only moderately changed. (An example is Gibbs sampling with small correlations.)
There exists a test mode that verifies whether the conditions for the method are satisfied or not. It can be switched on by calling unur_ars_set_verify and unur_ars_chg_verify respectively. Notice however that sampling is (much) slower then.
Method ARS aborts after a given number of iterations and return UNUR_INFINITY to prevent (almost) infinite loops. This might happen when the starting hat is much too large and it is not possible to insert new construction points due to severe numerical errors or (more likely) the given PDF is not log-concave. This maximum number of iterations can be set by means of a unur_ars_set_max_iter call.
Set maximum number of intervals. No construction points are added after the setup when the number of intervals suceeds max_ivs. It is increased automatically to twice the number of construction points if this is larger.
Default is
200.
Set construction points for the hat function. If cpoints is
NULLthen a heuristic rule of thumb is used to get n_cpoints construction points. This is the default behavior. n_cpoints should be at least2, otherwise defaults are used.The default number of construction points is 2.
— : int unur_ars_chg_reinit_percentiles (UNUR_GEN* generator, int n_percentiles, const double* percentiles)
By default, when the generator object is reinitialized, it used the same construction points as for the initialization procedure. Often the underlying distribution object has been changed only moderately. For example, the full conditional distribution of a multivariate distribution. In this case it might be more appropriate to use percentilesm of the hat function for the last (unchanged) distribution. percentiles must then be a pointer to an ordered array of numbers between
0.01and0.99. If percentiles isNULL, then a heuristic rule of thumb is used to get n_percentiles values for these percentiles. Notice that n_percentiles must be at least2, otherwise defaults are used. (Then the first and third quartiles are used by default.)
— : int unur_ars_chg_reinit_ncpoints (UNUR_GEN* generator, int ncpoints)
When reinit fails with the given construction points or the percentiles of the old hat function, another trial is undertaken with ncpoints construction points. ncpoints must be at least
10.Default:
30
The rejection loop stops after max_iter iterations and return UNUR_INFINITY.
Default:
10000
— : int unur_ars_chg_verify (UNUR_GEN* generator, int verify)
Turn verifying of algorithm while sampling on/off. If the condition squeeze(x) <= PDF(x) <= hat(x) is violated for some x then
unur_errnois set toUNUR_ERR_GEN_CONDITION. However notice that this might happen due to round-off errors for a few values of x (less than 1%).Default is
FALSE.
Sometimes it might happen that unur_init has been executed successfully. But when additional construction points are added by adaptive rejection sampling, the algorithm detects that the PDF is not log-concave.
With pedantic being
TRUE, the sampling routine is exchanged by a routine that simply returnsUNUR_INFINITY. Otherwise the new point is not added to the list of construction points. At least the hat function remains log-concave.Setting pedantic to
FALSEallows sampling from a distribution which is “almost” log-concave and small errors are tolerated. However it might happen that the hat function cannot be improved significantly. When the hat functions that has been constructed by the unur_init call is extremely large then it might happen that the generation times are extremely high (even hours are possible in extremely rare cases).Default is
FALSE.
Get the logarithm of area below the hat for the generator. (In case of an error
UNUR_INFINITYis returned.)
Evaluate the inverse of the CDF of the hat distribution at u.
If u is out of the domain [0,1] then
unur_errnois set toUNUR_ERR_DOMAINand the respective bound of the domain of the distribution are returned (which is-UNUR_INFINITYorUNUR_INFINITYin the case of unbounded domains).
Method CEXT is a wrapper for external generators for continuous univariate distributions. It allows the usage of external random variate generators within the UNU.RAN framework.
The following steps are required to use some external generator within the UNU.RAN framework (some of these are optional):
NULL. However, it is required if you want to pass
parameters of the generated distribution to the external
generator or for running some validation tests provided by
UNU.RAN.
int (*init)(UNUR_GEN *gen) and plug it into the generator
object using the
unur_cext_set_init
call. Notice that the
init routine must return UNUR_SUCCESS when it has
been executed successfully and UNUR_FAILURE otherwise.
It is possible to get the size of and the pointer to the array
of parameters of the underlying distribution object by the
respective calls
unur_cext_get_ndistrparams
and
unur_cext_get_distrparams.
Parameters for the external generator that are computed in the
init routine can be stored in a single array or structure
which is available by the
unur_cext_get_params
call.
Using an init routine is optional and can be omitted.
double (*sample)(UNUR_GEN *gen) and plug it into the
generator object using the
unur_cext_set_sample
call.
Uniform random numbers are provided by the unur_sample_urng call. Do not use your own implementation of a uniform random number generator directly. If you want to use your own random number generator we recommend to use the UNU.RAN interface (see see Using uniform random number generators).
The array or structure that contains parameters for the external generator that are computed in the init routine are available using the unur_cext_get_params call.
Using a sample routine is of course obligatory.
It is possible to change the parameters and the domain of the chosen distribution and run unur_reinit to reinitialize the generator object. The init routine is then called again.
Here is a short example that demonstrates the application of this method by means of the exponential distribution:
/* ------------------------------------------------------------- */
/* File: example_cext.c */
/* ------------------------------------------------------------- */
/* Include UNURAN header file. */
#include <unuran.h>
/* ------------------------------------------------------------- */
/* This example shows how an external generator for the */
/* exponential distribution with one scale parameter can be */
/* used within the UNURAN framework. */
/* */
/* Notice, that this example does not provide the simplest */
/* solution. */
/* ------------------------------------------------------------- */
/* Initialization routine. */
/* */
/* Here we simply read the scale parameter of the exponential */
/* distribution and store it in an array for parameters of */
/* the external generator. */
/* [ Of course we could do this in the sampling routine as */
/* and avoid the necessity of this initialization routine. ] */
int exponential_init (UNUR_GEN *gen)
{
/* Get pointer to parameters of exponential distribution */
double *params = unur_cext_get_distrparams(gen);
/* The scale parameter is the first entry (see manual) */
double lambda = (params) ? params[0] : 1.;
/* Get array to store this parameter for external generator */
double *genpar = unur_cext_get_params(gen, sizeof(double));
genpar[0] = lambda;
/* Executed successfully */
return UNUR_SUCCESS;
}
/* ------------------------------------------------------------- */
/* Sampling routine. */
/* */
/* Contains the code for the external generator. */
double exponential_sample (UNUR_GEN *gen)
{
/* Get scale parameter */
double *genpar = unur_cext_get_params(gen,0);
double lambda = genpar[0];
/* Sample a uniformly distributed random number */
double U = unur_sample_urng(gen);
/* Transform into exponentially distributed random variate */
return ( -log(1. - U) * lambda );
}
/* ------------------------------------------------------------- */
int main(void)
{
int i; /* loop variable */
double x; /* will hold the random number */
/* Declare the three UNURAN objects. */
UNUR_DISTR *distr; /* distribution object */
UNUR_PAR *par; /* parameter object */
UNUR_GEN *gen; /* generator object */
/* Use predefined exponential distribution with scale param. 2 */
double fpar[1] = { 2. };
distr = unur_distr_exponential(fpar, 1);
/* Use method CEXT */
par = unur_cext_new(distr);
/* Set initialization and sampling routines. */
unur_cext_set_init(par, exponential_init);
unur_cext_set_sample(par, exponential_sample);
/* Create the generator object. */
gen = unur_init(par);
/* It is important to check if the creation of the generator */
/* object was successful. Otherwise `gen' is the NULL pointer */
/* and would cause a segmentation fault if used for sampling. */
if (gen == NULL) {
fprintf(stderr, "ERROR: cannot create generator object\n");
exit (EXIT_FAILURE);
}
/* It is possible to reuse the distribution object to create */
/* another generator object. If you do not need it any more, */
/* it should be destroyed to free memory. */
unur_distr_free(distr);
/* Now you can use the generator object `gen' to sample from */
/* the standard Gaussian distribution. */
/* Eg.: */
for (i=0; i<10; i++) {
x = unur_sample_cont(gen);
printf("%f\n",x);
}
/* When you do not need the generator object any more, you */
/* can destroy it. */
unur_free(gen);
exit (EXIT_SUCCESS);
} /* end of main() */
/* ------------------------------------------------------------- */
Get default parameters for new generator.
Set initialization routine for external generator. Inside the
Important: The routine init must return
UNUR_SUCCESSwhen the generator was initialized successfully andUNUR_FAILUREotherwise.Parameters that are computed in the init routine can be stored in an array or structure that is avaiable by means of the unur_cext_get_params call. Parameters of the underlying distribution object can be obtained by the unur_cext_get_distrparams call.
Set sampling routine for external generator.
Important: Use
unur_sample_urng(gen)to get a uniform random number. The pointer to the array or structure that contains the parameters that are precomputed in the init routine are available byunur_cext_get_params(gen,0). Additionally one can use the unur_cext_get_distrparams call.
Get pointer to memory block for storing parameters of external generator. A memory block of size size is automatically (re-) allocated if necessary and the pointer to this block is stored in the generator object. If one only needs the pointer to this memory block set size to
0.Notice, that size is the size of the memory block and not the length of an array.
Important: This rountine should only be used in the initialization and sampling routine of the external generator.
— : int unur_cext_get_ndistrparams (UNUR_GEN* generator)
Get size of and pointer to array of parameters of underlying distribution in generator object.
Important: These rountines should only be used in the initialization and sampling routine of the external generator.
CSTD is a wrapper for special generators for continuous
univariate standard distributions. It only works for
distributions in the UNU.RAN library of standard distributions
(see Standard distributions).
If a distribution object is provided that is build from scratch,
it must provide the inverse CDF. Then CSTD implements the
inversion method. Otherwise, the NULL pointer is returned.
For some distributions more than one special generator is possible.
Create a distribution object for a standard distribution from the UNU.RAN library (see Standard distributions), or create a continuous distribution object and set the function for the inverse CDF using unur_distr_cont_set_invcdf. For some distributions more than one special generator (variants) is possible. These can be choosen by a unur_cstd_set_variant call. For possible variants see Standard distributions. However the following are common to all distributions:
UNUR_STDGEN_DEFAULTUNUR_STDGEN_FASTUNUR_STDGEN_INVERSIONNotice that the variant UNUR_STDGEN_FAST for a special
generator may be slower than one of the universal algorithms!
Additional variants may exist for particular distributions.
Sampling from truncated distributions (which can be constructed by changing the default domain of a distribution by means of unur_distr_cont_set_domain or unur_cstd_chg_truncated calls) is possible but requires the inversion method. Moreover the CDF of the distribution must be implemented.
It is possible to change the parameters and the domain of the chosen distribution and run unur_reinit to reinitialize the generator object.
Get default parameters for new generator. It requires a distribution object for a continuous univariant distribution from the UNU.RAN library of standard distributions (see Standard distributions).
Using a truncated distribution is allowed only if the inversion method is available and selected by the unur_cstd_set_variant call immediately after creating the parameter object. Use a unur_distr_cont_set_domain call to get a truncated distribution. To change the domain of a (truncated) distribution of a generator use the unur_cstd_chg_truncated call.
Set variant (special generator) for sampling from a given distribution. For possible variants see Standard distributions.
Common variants are
UNUR_STDGEN_DEFAULTfor the default generator,UNUR_STDGEN_FASTfor (one of the) fastest implemented special generators, andUNUR_STDGEN_INVERSIONfor the inversion method (if available). If the selected variant number is not implemented, then an error code is returned and the variant is not changed.
Change left and right border of the domain of the (truncated) distribution. This is only possible if the inversion method is used. Otherwise this call has no effect and an error code is returned.
Notice that the given truncated domain must be a subset of the domain of the given distribution. The generator always uses the intersection of the domain of the distribution and the truncated domain given by this call.
It is not required to run unur_reinit after this call has been used.
Important: If the CDF is (almost) the same for left and right and (almost) equal to
0or1, then the truncated domain is not chanced and the call returns an error code.Notice: If the parameters of the distribution has been changed it is recommended to set the truncated domain again, since the former call might change the domain of the distribution but not update the values for the boundaries of the truncated distribution.
HINV is a variant of numerical inversion, where the inverse CDF is approximated using Hermite interpolation, i.e., the interval [0,1] is split into several intervals and in each interval the inverse CDF is approximated by polynomials constructed by means of values of the CDF and PDF at interval boundaries. This makes it possible to improve the accuracy by splitting a particular interval without recomputations in unaffected intervals. Three types of splines are implemented: linear, cubic, and quintic interpolation. For linear interpolation only the CDF is required. Cubic interpolation also requires PDF and quintic interpolation PDF and its derivative.
These splines have to be computed in a setup step. However, it only works for distributions with bounded domain; for distributions with unbounded domain the tails are chopped off such that the probability for the tail regions is small compared to the given u-resolution.
The method is not exact, as it only produces random variates of the approximated distribution. Nevertheless, the maximal numerical error in "u-direction" (i.e. |U-CDF(X)|, for X = "approximate inverse CDF"(U) |U-CDF(X)|) can be set to the required resolution (within machine precision). Notice that very small values of the u-resolution are possible but may increase the cost for the setup step.
As the possible maximal error is only estimated in the setup it may be necessary to set some special design points for computing the Hermite interpolation to guarantee that the maximal u-error can not be bigger than desired. Such points are points where the density is not differentiable or has a local extremum. Notice that there is no necessity to do so. However, if you do not provide these points to the algorithm there might be a small chance that the approximation error is larger than the given u-resolution, or that the required number of intervals is larger than necessary.
HINV works for continuous univariate distribution objects with
given CDF and (optional) PDF. It uses Hermite interpolation of
order 1, 3 [default] or 5. The order can be set by means of
unur_hinv_set_order.
For distributions with unbounded domains the tails are chopped
off such that the probability for the tail regions is small
compared to the given u-resulution. For finding these cut points
the algorithm starts with the region [-1.e20,1.e20]. For
the exceptional case where this might be too small (or one knows
this region and wants to avoid this search heuristics) it can be
directly set via a
unur_hinv_set_boundary
call.
It is possible to use this method for generating from truncated distributions. It even can be changed for an existing generator object by an unur_hinv_chg_truncated call.
This method is not exact, as it only produces random variates of the approximated distribution. Nevertheless, the numerical error in "u-direction" (i.e. |U-CDF(X)|, for X = "approximate inverse CDF"(U) |U-CDF(X)|) can be controlled by means of unur_hinv_set_u_resolution. The possible maximal error is only estimated in the setup. Thus it might be necessary to set some special design points for computing the Hermite interpolation to guarantee that the maximal u-error can not be bigger than desired. Such points (e.g. extremal points of the PDF, points with infinite derivative) can be set using using the unur_hinv_set_cpoints call. If the mode for a unimodal distribution is set in the distribution object this mode is automatically used as design-point if the unur_hinv_set_cpoints call is not used.
As already mentioned the maximal error of this approximation is only estimated. If this error is crucial for an application we recommend to compute this error using unur_hinv_estimate_error which runs a small Monte Carlo simulation.
It is possible to change the parameters and the domain of the chosen distribution and run unur_reinit to reinitialize the generator object. The values given by the last unur_hinv_chg_truncated call will be then changed to the values of the domain of the underlying distribution object. Moreover, starting construction points (nodes) that are given by a unur_hinv_set_cpoints call are ignored when unur_reinit is called. It is important to note that for a distribution from the UNU.RAN library of standard distributions (see Standard distributions) the normalization constant has to be updated using the unur_distr_cont_upd_pdfarea call whenever its parameters have been changed by means of a unur_distr_cont_set_pdfparams call.
Set order of Hermite interpolation. Valid orders are
1,3, and5. Notice that order greater than1requires the density of the distribution, and order greater than3even requires the derivative of the density. Using order1results for most distributions in a huge number of intervals and is therefore not recommended. If the maximal error in u-direction is very small (say smaller than1.e-10), order5is recommended as it leads to considerably fewer design points, as long there are no poles or heavy tails.Remark: When the target distribution has poles or (very) heavy tails order
5(i.e., quintic interpolation) is numerically less stable and more sensitive to round-off errors than order3(i.e., cubic interpolation).Default is
3if the density is given and1otherwise.
Set maximal error in u-direction. However, the given u-error must not be smaller than machine epsilon (
DBL_EPSILON) and should not be too close to this value. As the resolution of most uniform random number sources is 2^(-32) =2.3e-10, a value of1.e-10leads to an inversion algorithm that could be called exact. For most simulations slightly bigger values for the maximal error are enough as well.Remark: The u-error might become larger than u_resolution due to rescaling of floating point numbers when the domain of the distribution is truncated by a unur_hinv_chg_truncated call.
Default is
1.e-10.
Set starting construction points (nodes) for Hermite interpolation.
As the possible maximal error is only estimated in the setup it may be necessary to set some special design points for computing the Hermite interpolation to guarantee that the maximal u-error can not be bigger than desired. We suggest to include as special design points all local extrema of the density, all points where the density is not differentiable, and isolated points inside of the domain with density 0. If there is an interval with density constant equal to 0 inside of the given domain of the density, both endpoints of this interval should be included as special design points. Notice that there is no necessity to do so. However, if these points are not provided to the algorithm the approximation error might be larger than the given u-resolution, or the required number of intervals could be larger than necessary.
Important: Notice that the given points must be in increasing order and they must be disjoint.
Important: The boundary point of the computational region must not be given in this list! Points outside the boundary of the computational region are ignored.
Default is for unimodal densities - if known - the mode of the density, if it is not equal to the border of the domain.
Set the left and right boundary of the computational interval. Of course
+/- UNUR_INFINITYis not allowed. If the CDF at left and right is not close to the respective values0.and1.then this interval is increased by a (rather slow) search algorithm.Important: This call does not change the domain of the given distribution itself. But it restricts the domain for the resulting random variates.
Default is
1.e20.
Set factor for relative size of the guide table for indexed search (see also method DGT DGT). It must be greater than or equal to
0. When set to0, then sequential search is used.Default is
1.
Set maximum number of intervals. No generator object is created if the necessary number of intervals for the Hermite interpolation exceeds max_ivs. It is used to prevent the algorithm to eat up all memory for very badly shaped CDFs.
Default is
1000000(1.e6).
Get number of nodes (design points) used for Hermite interpolation in the generator object. The number of intervals is the number of nodes minus 1. It returns an error code in case of an error.
Evaluate Hermite interpolation of inverse CDF at u. If u is out of the domain [0,1] then
unur_errnois set toUNUR_ERR_DOMAINand the respective bound of the domain of the distribution are returned (which is-UNUR_INFINITYorUNUR_INFINITYin the case of unbounded domains).Notice: When the domain has been truncated by a unur_hinv_chg_truncated call then the inverse CDF of the truncated distribution is returned.
Changes the borders of the domain of the (truncated) distribution.
Notice that the given truncated domain must be a subset of the domain of the given distribution. The generator always uses the intersection of the domain of the distribution and the truncated domain given by this call. The tables of splines are not recomputed. Thus it might happen that the relative error for the generated variates from the truncated distribution is greater than the bound for the non-truncated distribution. This call also fails when the CDF values of the boundary points are too close, i.e. when only a few different floating point numbers would be computed due to round-off errors with floating point arithmetic.
Remark: The u-error might become larger than the u_resolution given by a unur_hinv_set_u_resolution call due to rescaling of floating point numbers when the domain of the distribution is truncated.
When failed an error code is returned.
Important: Always check the return code since the domain is not changed in case of an error.
Estimate maximal u-error and mean absolute error (MAE) for generator by means of a (quasi-) Monte-Carlo simulation with sample size samplesize. The results are stored in max_error and MAE, respectively.
It returns
UNUR_SUCCESSif successful.
Generates random variate with given hazard rate which must be bounded from above. It uses the thinning method with a constant dominating hazard function.
HRB requires a hazard function for a continuous distribution together with an upper bound. The latter has to be set using the unur_hrb_set_upperbound call. If no such upper bound is given it is assumed that the upper bound can be achieved by evaluating the hazard rate at the left hand boundary of the domain of the distribution. Notice, however, that for decreasing hazard rate the method HRD (see Hazard Rate Decreasing) is much faster and thus the prefered method.
It is important to note that the domain of the distribution can
be set via a
unur_distr_cont_set_domain
call.
However, the left border must not be negative. Otherwise it is
set to 0. This is also the default if no domain is
given at all. For computational reasons the right border is
always set to UNUR_INFINITY independently of the given
domain. Thus for domains bounded from right the function for
computing the hazard rate should return UNUR_INFINITY
right of this domain.
For distributions with increasing hazard rate method HRI (see Hazard Rate Increasing) is required.
It is possible to change the parameters and the domain of the chosen distribution and run unur_reinit to reinitialize the generator object. Notice, that the upper bound given by the unur_hrb_set_upperbound call cannot be changed and must be valid for the changed distribution.
Set upper bound for hazard rate. If this call is not used it is assumed that the the maximum of the hazard rate is achieved at the left hand boundary of the domain of the distribution.
— : int unur_hrb_chg_verify (UNUR_GEN* generator, int verify)
Turn verifying of algorithm while sampling on/off. If the hazard rate is not bounded by the given bound, then
unur_errnois set toUNUR_ERR_GEN_CONDITION.Default is
FALSE.
Generates random variate with given non-increasing hazard rate. It is necessary that the distribution object contains this hazard rate. Decreasing hazard rate implies that the corresponding PDF of the distribution has heavier tails than the exponential distribution (which has constant hazard rate).
HRD requires a hazard function for a continuous distribution with non-increasing hazard rate. There are no parameters for this method.
It is important to note that the domain of the distribution can
be set via a
unur_distr_cont_set_domain
call. However, only
the left hand boundary is used. For computational reasons the
right hand boundary is always reset to UNUR_INFINITY.
If no domain is given by the user then the left hand boundary is
set to 0.
For distributions which do not have decreasing hazard rates but are bounded from above use method HRB (see Hazard Rate Bounded). For distributions with increasing hazard rate method HRI (see Hazard Rate Increasing) is required.
It is possible to change the parameters and the domain of the chosen distribution and run unur_reinit to reinitialize the generator object.
— : int unur_hrd_chg_verify (UNUR_GEN* generator, int verify)
Turn verifying of algorithm while sampling on/off. If the hazard rate is not bounded by the given bound, then
unur_errnois set toUNUR_ERR_GEN_CONDITION.Default is
FALSE.
Generates random variate with given non-increasing hazard rate. It is necessary that the distribution object contains this hazard rate. Increasing hazard rate implies that the corresponding PDF of the distribution has heavier tails than the exponential distribution (which has constant hazard rate).
The method uses a decomposition of the hazard rate into a main part which is constant for all x beyond some point p0 and a remaining part. From both of these parts points are sampled using the thinning method and the minimum of both is returned. Sampling from the first part is easier as we have a constant dominating hazard rate. Thus p0 should be large. On the other hand, if p0 is large than the thinning algorithm needs many iteration. Thus the performance of the the algorithm deponds on the choice of p0. We found that values close to the expectation of the generated distribution result in good performance.
HRI requires a hazard function for a continuous distribution with non-decreasing hazard rate. The parameter p0 should be set to a value close to the expectation of the required distribution using unur_hri_set_p0. If performance is crucial one may try other values as well.
It is important to note that the domain of the distribution can
be set via a
unur_distr_cont_set_domain
call. However, only
the left hand boundary is used. For computational reasons the
right hand boundary is always reset to UNUR_INFINITY.
If no domain is given by the user then the left hand boundary is
set to 0.
For distributions with decreasing hazard rate method HRD (see Hazard Rate Decreasing) is required. For distributions which do not have increasing or decreasing hazard rates but are bounded from above use method HRB (see Hazard Rate Bounded).
It is possible to change the parameters and the domain of the chosen distribution and run unur_reinit to reinitialize the generator object.
Notice, that the upper bound given by the unur_hrb_set_upperbound call cannot be changed and must be valid for the changed distribution. Notice that the parameter p0 which has been set by a unur_hri_set_p0 call cannot be changed and must be valid for the changed distribution.
Set design point for algorithm. It is used to split the domain of the distribution. Values for p0 close to the expectation of the distribution results in a relatively good performance of the algorithm. It is important that the hazard rate at this point must be greater than
0and less thanUNUR_INFINITY.Default: left boundary of domain +
1.
— : int unur_hri_chg_verify (UNUR_GEN* generator, int verify)
Turn verifying of algorithm while sampling on/off. If the hazard rate is not bounded by the given bound, then
unur_errnois set toUNUR_ERR_GEN_CONDITION.Default is
FALSE.
ITDR is an acceptance/rejection method that works for monotone densities. It is especially designed for PDFs with a single pole. It uses different hat functions for the pole region and for the tail region. For the tail region Transformed Density Rejection with a single construction point is used. For the pole region a variant called Inverse Transformed Density Rejection is used. The optimal splitting point between the two regions and the respective maximum local concavity and inverse local concavity (see Glossary) that guarantee valid hat functions for each regions are estimated. This splitting point is set to the intersection point of local concavity and inverse local concavity. However, it is assumed that both, the local concavity and the inverse local concavity do not have a local minimum in the interior of the domain (which is the case for all standard distributions with a single pole). In other cases (or when the built-in search routines do not compute non-optimal values) one can provide the splitting point, and the c-values.
Method ITDR requires a distribution object with given PDF and its derivative and the location of the pole (or mode). The PDF must be monotone and may contain a pole. It must be set via the unur_distr_cont_set_pdf and unur_distr_cont_set_dpdf calls. The PDF should return UNUR_INFINITY for the pole. Alternatively, one can also set the logarithm of the PDF and its derivative via the unur_distr_cont_set_logpdf and unur_distr_cont_set_dlogpdf calls. This is in especially useful since then the setup and search routines are numerically more stable. Moreover, for many distributions computing the logarithm of the PDF is less expensive then computing the PDF directly.
The pole of the distribution is given by a
unur_distr_cont_set_mode
call. Notice that distributions with
“heavy” poles may have numerical problems caused by the
resultion of the floating point numbers used by computers.
While the minimal distance between two different floating point
numbers is about 1.e-320 near 0. it increases
to 1.e-16 near 1. Thus any random variate
generator implemented on a digital computer in fact draws samples
from a discrete distribution that approximates the desired
continuous distribution. For distributions with “heavy” poles
not at 0 this approximation may be too crude and thus every
goodness-of-fit test will fail.
Besides this theoretic problem that cannot be resolved we
have to take into consideration that round-off errors occur more
frequently when we have PDFs with poles far away from
0. Method ITDR tries to handles this situation as good as
possible by moving the pole into 0.
Thus do not use a wrapper for your PDF that hides this shift
since the information about the resolution of the floating point
numbers near the pole gets lost.
Method ITDR uses different hats for the pole region and for the
tail region. The splitting point between these two regions, the
optimal c-value and design points for constructing the hats
using Transformed Density Rejection are computed automatically.
(The results of these computations can be read using the
respective calls
unur_itdr_get_xi
unur_itdr_get_cp
, and
unur_itdr_get_ct
for the intersection point between local
concavity and inverse local concavity, the c-value for the
pole and the tail region.)
However, one can also analyze the local concavity and inverse
local concavity set the corresponding values using
unur_itdr_set_xi
unur_itdr_set_cp
, and
unur_itdr_set_ct
calls.
Notice, that c-values greater than -1/2 can be set to
-0.5. Although this results in smaller acceptance
probabities sampling from the hat distribution is much faster
than for other values of c. Depending on the expenses of
evaluating the PDF the resulting algorithm is usually faster.
It is possible to change the parameters and the domain of the chosen distribution and run unur_reinit to reinitialize the generator object. However, the values given by unur_itdr_set_xi unur_itdr_set_cp , or unur_itdr_set_ct calls are then ignored when unur_reinit is called.
Sets points where local concavity and inverse local concavity are (almost) equal. It is used to estimate the respective c-values for pole region and hat regions and to determine the splitting point bx between pole and tail region. If no such point is provided it will be computed automatically.
Default: not set.
Sets parameter cp for transformation T for inverse density in pole region. It must be at most 0 and greater than -1. A value of
-0.5is treated separately and usually results in faster marginal generation time (at the expense of smaller acceptance probabilities. If no cp-value is given it is estimated automatically.Default: not set.
Sets parameter ct for transformation T for density in tail region. It must be at most 0. For densities with unbounded domain it must be greater than -1. A value of
-0.5is treated separately and usually results in faster marginal generation time (at the expense of smaller acceptance probabilities. If no ct-value is given it is estimated automatically.Default: not set.
— : double unur_itdr_get_cp (UNUR_GEN* generator)
— : double unur_itdr_get_ct (UNUR_GEN* generator)
Get intersection point xi, and c-values cp and ct, respectively. (In case of an error
UNUR_INFINITYis returned.)
Get area below hat. (In case of an error
UNUR_INFINITYis returned.)
Turn verifying of algorithm while sampling on/off.
If the condition PDF(x) <= hat(x)is violated for some x then
unur_errnois set toUNUR_ERR_GEN_CONDITION. However, notice that this might happen due to round-off errors for a few values of x (less than 1%).Default is
FALSE.
Change the verifying of algorithm while sampling on/off.
NINV implementations of some methods for numerical inversion: Newton's method, regula falsi (combined with interval bisectioning), and bisection method. Regula falsi and bisection method require only the CDF while Newton's method also requires the PDF. To speed up marginal generation times a table with suitable starting points can be created during the setup. The performance of the algorithm can adjusted by the desired accuracy of the method. It is possible to use this method for generating from truncated distributions. The truncated domain can be changed for an existing generator object.
Method NINV generates random variates by numerical inversion and requires a continuous univariate distribution objects with given CDF. Three variants are available:
Newton's method additionally requires the PDF of the distribution and cannot be used otherwise (NINV automatically switches to regula falsi then). Default algorithm is regula falsi. It is slightly slower but numerically much more stable than Newton's algorithm. Interval bisectioning is the slowest method and should only be considered as a last resort when the other methods fails.
It is possible to draw samples from truncated distributions. The truncated domain can even be changed for an existing generator object by an unur_ninv_chg_truncated call.
Marginal generation times can be sped up by means of a table with suitable starting points which can be created during the setup. Using such a table can be switched on by means of a unur_ninv_set_table call where the table size is given as a parameter. The table is still useful when the (truncated) domain is changed often, since it is computed for the domain of the given distribution. (It is not possible to enlarge this domain.) If it is necessary to recalculate the table during sampling, the command unur_ninv_chg_table can be used. As a rule of thumb using such a table is appropriate when the number of generated points exceeds the table size by a factor of 100.
The default number of iterations of NINV should be enough for all reasonable cases. Nevertheless, it is possible to adjust the maximal number of iterations with the commands unur_ninv_set_max_iter and unur_ninv_chg_max_iter. In particular this might be necessary when the PDF has a pole or the distribution has extremely heavy tails.
It is also possible to set/change the accuracy of the method (which also heavily influencies the generation time). We use two measures for the approximation error which can be used independently: x-error and u-error (see Inversion for more details). It is possible to set the maximal tolerated error using with unur_ninv_set_x_resolution and with unur_ninv_set_u_resolution resp., and change it with the respective calls unur_ninv_chg_x_resolution and unur_ninv_chg_x_resolution. The algorithm tries to satisfy both accuracy goals (and raises an error flag it this fails). One of these accuracy checks can be disabled by setting the accuracy goal to a negative value.
NINV tries to use proper starting values for both the regula falsi and bisection method, and for Newton's method. Of course the user might have more knowledge about the properties of the target distribution and is able to share his wisdom with NINV using the respective commands unur_ninv_set_start and unur_ninv_chg_start. It is possible to change the parameters and the domain of the chosen distribution and run unur_reinit to reinitialize the generator object. The values given by the last unur_ninv_chg_truncated call will be then changed to the values of the domain of the underlying distribution object. It is important to note that for a distribution from the UNU.RAN library of standard distributions (see Standard distributions) the normalization constant has to be updated using the unur_distr_cont_upd_pdfarea call whenever its parameters have been changed by means of a unur_distr_cont_set_pdfparams call.
It might happen that NINV aborts
unur_sample_cont
without
computing the correct value (because the maximal number
iterations has been exceeded). Then the last approximate value
for x is returned (with might be fairly false) and
unur_error is set to UNUR_ERR_GEN_SAMPLING.
Switch to regula falsi combined with interval bisectioning. (This the default.)
Switch to Newton's method. Notice that it is numerically less stable than regula falsi. It it is not possible to invert the CDF for a particular uniform random number U when calling unur_sample_cont
unur_erroris set toUNUR_ERR_GEN_SAMPLING. Thus it is recommended to checkunur_errorbefore using the result of the sampling routine.
Switch to bisection method. This is a slow algorithm and should only be used as a last resort.
— : int unur_ninv_chg_max_iter (UNUR_GEN* generator, int max_iter)
Set and change number of maximal iterations. Default is
100.
— : int unur_ninv_chg_x_resolution (UNUR_GEN* generator, double x_resolution)
Set and change the maximal tolerated relative x-error. If x_resolution is negative then checking of the x-error is disabled.
Default is
1.e-8.
— : int unur_ninv_chg_u_resolution (UNUR_GEN* generator, double u_resolution)
Set and change the maximal tolerated (abolute) u-error. If u_resolution is negative then checking of the u-error is disabled.
Default is
-1(disabled).
Set starting points. If not set, suitable values are chosen automatically.
Newton: left: starting point Regula falsi: left, right: boundary of starting interval If the starting points are not set then the follwing points are used by default:
Newton: left: CDF(left) = 0.5 Regula falsi: left: CDF(left) = 0.1 right: CDF(right) = 0.9 If left == right, then UNU.RAN always uses the default starting points!
Change the starting points for numerical inversion. If left==right, then UNU.RAN uses the default starting points (see unur_ninv_set_start ).
Generates a table with no_of_points points containing suitable starting values for the iteration. The value of no_of_points must be at least 10 (otherwise it will be set to 10 automatically).
The table points are chosen such that the CDF at these points form an equidistance sequence in the interval (0,1).
If a table is used, then the starting points given by unur_ninv_set_start are ignored.
No table is used by default.
Recomputes a table as described in unur_ninv_set_table.
Changes the borders of the domain of the (truncated) distribution.
Notice that the given truncated domain must be a subset of the domain of the given distribution. The generator always uses the intersection of the domain of the distribution and the truncated domain given by this call. Moreover the starting point(s) will not be changed.
Important: If the CDF is (almost) the same for left and right and (almost) equal to
0or1, then the truncated domain is not chanced and the call returns an error code.Notice: If the parameters of the distribution has been changed by a unur_distr_cont_set_pdfparams call it is recommended to set the truncated domain again, since the former call might change the domain of the distribution but not update the values for the boundaries of the truncated distribution.
Get approximate approximate value of inverse CDF at u. If u is out of the domain [0,1] then
unur_errnois set toUNUR_ERR_DOMAINand the respective bound of the domain of the distribution are returned (which is-UNUR_INFINITYorUNUR_INFINITYin the case of unbounded domains).Notice: This function always evaluates the inverse CDF of the given distribution. A call to unur_ninv_chg_truncated call has no effect.
NROU is an implementation of the (generalized) ratio-of-uniforms method which uses (minimal) bounding rectangles, see Ratio-of-Uniforms. It uses a positive control parameter r for adjusting the algorithm to the given distribution to improve performance and/or to make this method applicable. Larger values of r increase the class of distributions for which the method works at the expense of a higher rejection constant. For computational reasons r=1 should be used if possible (this is the default). Moreover, this implementation uses the center muof the distribution (see unur_distr_cont_get_center for details of its default values).
For the special case with r=1the coordinates of the minimal bounding rectangles are given by
v+ = sup_x sqrt(PDF(x)),
u- = inf_x (x- mu) sqrt(PDF(x)),
u+ = sup_x (x- mu) sqrt(PDF(x)),
where muis the center of the distribution. For other values of r we have
v+ = sup_x (PDF(x))1/(r+1),
u- = inf_x (x- mu) (PDF(x))r/(r+1),
u+ = sup_x (x- mu) (PDF(x))r/(r+1).
These bounds can be given directly. Otherwise they are computed automatically by means of a (slow) numerical routine. Of course this routine can fail, especially when this rectangle is not bounded.
It is important to note that the algorithm works with PDF(x- mu)instead of PDF(x). This is important as otherwise the acceptance region can become a very long and skinny ellipsoid along a diagonal of the (huge) bounding rectangle.
For using the NROU method UNU.RAN needs the PDF of the distribution. Additionally, the parameter r can be set via a unur_vnrou_set_r call. Notice that the acceptance probability decreases when r is increased. On the other hand is is more unlikely that the bounding rectangle does not exist if r is small.
A bounding rectangle can be given by the unur_vnrou_set_u and unur_vnrou_set_v calls.
Important: The bounding rectangle has to be
provided for the function
PDF(x-center)! Notice that center is the center of the given
distribution, see
unur_distr_cont_set_center.
If in doubt or if this value is not optimal, it can be changed
(overridden) by a
unur_nrou_set_center
call.
If the coordinates of the bounding rectangle are not provided by the user then the minimal bounding rectangle is computed automatically.
By means of unur_vnrou_set_verify and unur_vnrou_chg_verify one can run the sampling algorithm in a checking mode, i.e., in every cycle of the rejection loop it is checked whether the used rectangle indeed enclosed the acceptance region of the distribution. When in doubt (e.g., when it is not clear whether the numerical routine has worked correctly) this can be used to run a small Monte Carlo study.
It is possible to change the parameters and the domain of the chosen distribution and run unur_reinit to reinitialize the generator object. Notice, that derived parameters like the mode must also be (re-) set if the parameters or the domain has be changed. Notice, however, that then the values that has been set by unur_vnrou_set_u and unur_vnrou_set_v calls are removed and the coordinates of the bounding box are computed numerically.
Sets left and right boundary of bounding rectangle. If no values are given, the boundary of the minimal bounding rectangle is computed numerically.
Notice: Computing the minimal bounding rectangle may fail under some circumstances. Moreover, for multimodal distributions the bounds might be too small as only local extrema are computed. Nevertheless, for T_c-concave distributions with c=-1/2it should work.
Important: The bounding rectangle that has to be provided is for the function PDF(x-center)!
Default: not set.
Set upper boundary for bounding rectangle. If this value is not given then sqrt(PDF(mode))is used instead.
Notice: When the mode is not given for the distribution object, then it will be computed numerically.
Default: not set.
Sets the parameter r of the generalized ratio-of-uniforms method.
Notice: This parameter must satisfy r>0.
Default:
1.
Set the center (mu) of the PDF. If not set the center of the given distribution object is used.
Default: see unur_distr_cont_set_center.
Turn verifying of algorithm while sampling on/off.
If the condition PDF(x) <= hat(x)is violated for some x then
unur_errnois set toUNUR_ERR_GEN_CONDITION. However, notice that this might happen due to round-off errors for a few values of x (less than 1%).Default is
FALSE.
Change the verifying of algorithm while sampling on/off.
PINV is a variant of numerical inversion, where the inverse CDF is approximated using Newton's interpolating formula. The interval [0,1] is split into several subintervals. In each of these the inverse CDF is constructed at nodes (CDF(x),x)for some points x in this subinterval. If the PDF is given then the CDF is computed numerically from the given PDF using adaptive Gauss-Lobatto integration with 5 points. Subintervals are splitted until the requested accuracy goal is reached.
The method is not exact, as it only produces random variates of the approximated distribution. Nevertheless, the maximal tolerated approximation error can be set to be the resolution (but of course is bounded by the machine precision). We use the u-error |U-CDF(X)|to measure the error for X = "approximate inverse CDF"(U). Notice that very small values of the u-resolution are possible but increase the cost for the setup step. We call the maximal tolerated u-error the u-resolution of the algorithm in the sequel.
Both the order of the interpolating polynomial and the u-resolution can be selected.
The interpolating polynomials have to be computed in a setup step. However, it only works for distributions with bounded domain; for distributions with unbounded domain the tails are cut off such that the probability for the tail regions is small compared to the given u-resolution.
The construction of the interpolation polynomial only works when the PDF is unimodal or when the PDF does not vanish between two modes.
There are some restrictions for the given distribution:
1.e-12.
PINV works for continuous univariate distribution objects with
given PDF. The corresponding distribution object must contain a
typical point of the distribution, i.e., a point where the PDF
is not too small, e.g., (a point near) the mode.
However, it is important that the center is not the
pole of the distribution (or a point too close to the pole).
It can be set using a
unur_distr_cont_set_center
or
a
unur_distr_cont_set_mode
call. (If neither is set, 0
is assumed!)
It is recommended that the domain of the distribution with
bounded domain is specified using a
unur_distr_cont_set_domain
call. Otherwise, the boundary is searched numerically which
might be rather expensive, especially when this boundary point
is 0.
When sampling from truncated distributions with extreme truncation points, it is recommended to provide the log-density using unur_distr_cont_set_logpdf and the mode. Then the PDF is rescaled such that the PDF at the mode is 1. Thus the algorithm is numerically more stable.
It is also possible to use the CDF of the distribution instead of the PDF. Then the distribution object must contain a pointer to the CDF. Moreover, this variant of the algorithmus has to be switched on using an unur_pinv_set_usecdf call. Notice, however, that the setup for this variant is numerically less stable than using integration of the PDF (the default variant).
The inverse CDF is interpolated using Newton polymials. The order of this polynomial can be set by means of a unur_pinv_set_order call.
For distributions with unbounded domains the tails are cut
off such that the probability for the tail regions is small
compared to the given u-resulution. For finding these cut points
the algorithm starts with the region [-1.e100,1.e100]. For
the exceptional case where this does not work these starting
points can be changed via a
unur_pinv_set_boundary
call.
This method is not exact, as it only produces random variates of the approximated distribution. Nevertheless, the numerical error in "u-direction" (i.e., |U-CDF(X)|, for X = "approximate inverse CDF"(U) |U-CDF(X)|) can be controlled by means of unur_pinv_set_u_resolution. However, the maximal error of this approximation is only estimated. For very small u-resolutions the actual approximation error might be (slightly) larger than the requested u-resolution. (Of course the size of this value depends on the given PDF.) If this error is crucial for an application we recommend to compute this error using unur_pinv_estimate_error which runs a small Monte Carlo simulation. See also the documentation for function unur_pinv_set_u_resolution and the remark given there.
The number of required subintervals heavily depends on the order of the interpolating polynomial and the requested u-resolution: it increases when order or u-resolution are decreased. It can be checked using a unur_pinv_get_n_intervals call. The maximum number of such subintervals is fixed but can be increased using a unur_pinv_set_max_intervals call. If this maximum number is too small then the set-up aborts with a corresponding error message.
Set order of interpolation. Valid orders are between
3and12. Higher orders result in fewer intervals for the approximations.Default:
5.
Set maximal tolerated u-error. Values of u_resolution must at least
1.e-15and1.e-5at most. Notice that the resolution of most uniform random number sources is 2-32=2.3e-10. Thus a value of1.e-10leads to an inversion algorithm that could be called exact. For most simulations slightly bigger values for the maximal error are enough as well.Smaller values for u_resolution increase the number of subinterval that are necessary for the approximation of the inverse CDF. For very small values (less then
1.e-12) this number might exceed the maximum number of such intervals. However, this number can be increased using a unur_pinv_set_max_intervals call.Remark: We ran many experiments and found that the observed u-error was always smaller than the given u_resolution whenever this value was
1.e-12. For values smaller than1e-13the maximal observed u-error was slightly larger. One use1.e-15if best approximation is required. However, then the actual u-error can be as large as1.e-14.Warning! These figures are based on our experiments (with some tolarence added to be on the safe side). There is no guarentee for these error estimates for a particular distribution.
Default is
1.e-10.
Use PDF (if available) to compute approximate inverse CDF.
This is the default.
Use CDF (if available) to compute approximate inverse CDF.
Remark: We ran many experiments and found that for small values of the given u_resolution (less than
1.e-12) the setup fails for distributions with heavy tails. We found that using the PDF (instead of the CDF) is numerically more stable.
Set left and right point for finding the cut-off points for the "computational domain", i.e., the domain that covers the essential part of the distribution. The cut-off points are computed such that the tail probabilities are smaller than given by unur_pinv_set_u_resolution. It is usually safe to use a large interval. However,
+/- UNUR_INFINITYis not allowed.Important: This call does not change the domain of the given distribution itself. But it restricts the domain for the resulting random variates.
Default: intersection of
[-1.e100,+1.e100]and the given domain of the distribution.
If left or right is set to
FALSEthen the respective boundary as given by a unur_pinv_set_boundary call is used without any further computations. However, these boundary points might cause numerical problems during the setup when PDF returns0“almost everywhere”. If set toTRUE(the default) then the computational interval is shortened to a more sensible region by means of a search algorithm. Switching off this search is useful, e.g., for the Gamma(2) distribution where the left border0is fixed and finite.Remark: The searching algorithm assumes that the support of the distribution is connected.
Remark: Do not set this parameter to
FALSEexcept when searching for cut-off points fails and one wants to try with precomputed values.Default:
TRUE.
Set maximum number of intervals. max_ivs must be at least
100and at most1000000.Default is
10000.
Get number of intervals used for interpolation in the generator object. It returns
0in case of an error.
Evaluate interpolation of inverse CDF at u. If u is out of the domain (0,1) then
unur_errnois set toUNUR_ERR_DOMAINand the respective bound of the domain of the distribution are returned (which is-UNUR_INFINITYorUNUR_INFINITYin the case of unbounded domains).
Estimate maximal u-error and mean absolute error (MAE) for generator by means of Monte-Carlo simulation with sample size samplesize. The results are stored in max_error and MAE, respectively.
It returns
UNUR_SUCCESSif successful.
SROU is based on the ratio-of-uniforms method (see Ratio-of-Uniforms) that uses universal inequalities for constructing a (universal) bounding rectangle. It works for all T-concave distributions, including log-concave and T-concave distributions with T(x) = -1/sqrt(x).
Moreover an (optional) parameter r can be given, to
adjust the generator to the given distribution. This parameter
is strongly related to the parameter c for transformed
density rejection (see TDR) via the formula
c = -r/(r+1). The rejection constant increases with higher
values for r. On the other hand, the given density must
be
T_c-concave for the corresponding c.
The default setting for r is 1 which results in a very
simple code. (For other settings, sampling uniformly from the
acceptance region is more complicated.)
Optionally the CDF at the mode can be given to increase the
performance of the algorithm. Then the rejection constant is
reduced by 1/2 and (if r=1) even a universal squeeze can
(but need not be) used.
A way to increase the performance of the algorithm when the
CDF at the mode is not provided is the usage of the mirror
principle (only if r=1). However, using squeezes and using
the mirror principle is only recommended when the PDF is
expensive to compute.
The exact location of the mode and/or the area below the PDF can be replace by appropriate bounds. Then the algorithm still works but has larger rejection constants.
SROU works for any continuous univariate distribution object with
given
T_c-concave PDF with
c<1,)
mode and area below PDF. Optional the CDF at the mode
can be given to increase the performance of the algorithm by
means of the
unur_srou_set_cdfatmode
call. Additionally
squeezes can be used and switched on via
unur_srou_set_usesqueeze
(only if r=1).
A way to increase the performance of the algorithm when the
CDF at the mode is not provided is the usage of the mirror
principle which can be swithced on by means of a
unur_srou_set_usemirror
call (only if r=1) .
However using squeezes and using
the mirror principle is only recommended when the PDF is
expensive to compute.
The parameter r can be given, to adjust the generator to
the given distribution. This parameter is strongly related
parameter c for transformed density rejection via the
formula c = -r/(r+1).
The parameter r can be any value larger than or equal to
1. Values less then 1 are automatically set to 1.
The rejection constant depends on the chosen parameter
r but not on the particular distribution. It is 4 for
r equal to 1 and higher for higher values of r.
It is important to note that different algorithms for different
values of r: If r equal to 1 this is much faster
than the algorithm for r greater than 1.
The default setting for r is 1.
If the (exact) area below the PDF is not known, then an upper bound can be used instead (which of course increases the rejection constant). But then the squeeze flag must not be set and unur_srou_set_cdfatmode must not be used.
If the exact location of the mode is not known, then use the approximate location and provide the (exact) value of the PDF at the mode by means of the unur_srou_set_pdfatmode call. But then unur_srou_set_cdfatmode must not be used. Notice, that a (slow) numerical mode finder will be used if no mode is given at all. It is even possible to give an upper bound for the PDF only. However, then the (upper bound for the) area below the PDF has to be multiplied by the ratio between the upper bound and the lower bound of the PDF at the mode. Again setting the squeeze flag and using unur_srou_set_cdfatmode is not allowed.
It is possible to change the parameters and the domain of the chosen distribution and run unur_reinit to reinitialize the generator object. Notice, that derived parameters like the mode must also be (re-) set if the parameters or the domain has be changed. Moreover, if the PDF at the mode has been provided by a unur_srou_set_pdfatmode call, additionally unur_srou_chg_pdfatmode must be used (otherwise this call is not necessary since then this figure is computed directly from the PDF).
There exists a test mode that verifies whether the conditions for the method are satisfied or not while sampling. It can be switched on by calling unur_srou_set_verify and unur_srou_chg_verify respectively. Notice however that sampling is (a little bit) slower then.
Set parameter r for transformation. Only values greater than or equal to 1 are allowed. The performance of the generator decreases when r is increased. On the other hand r must not be set to small, since the given density must be T_c-concave for c = -r/(r+1).
Notice: If r is set to
1a simpler and much faster algorithm is used then for r greater than one.For computational reasons values of r that are greater than
1but less than1.01are always set to1.01.Default is
1.