pseudorandom.h

随机数类,产生随机数

/* All the sampling source code are taken from GSL - GNU Scientific Library 
 * except multigaussian distribution.i copy it from  
 * John Burkardt's  PROB routines with my modification
 * 
 * Copyright (C) 1996, 1997, 1998, 1999, 2000 James Theiler, Brian Gough
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or (at
 * your option) any later version.
 * 
 * This program is distributed in the hope that it will be useful, but
 * WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * General Public License for more details.
 * 
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
 */
# ifndef PseudoRandom_H_
# define PseudoRandom_H_

# include "MersenneTwister.h"
# include "matrix.hpp"

# ifndef M_E
# define M_E        2.71828182845904523536028747135      /* e */
# endif

# ifndef M_PI
# define M_PI       3.14159265358979323846264338328      /* pi */
# endif

using namespace ldas;

// C++ interface to the random-distribution generators
class PseudoRandom : public MTRand
{
 public:
	PseudoRandom() : MTRand() {};  // Auto-initialize with /dev/urandom or time() and clock()
	PseudoRandom(unsigned long seed) : MTRand( seed ) {};  // Initialize with a given seed
	PseudoRandom(unsigned long * array, unsigned long length = MTRand::N) : MTRand( array, length ) {};  // Initialize by an array with array-length
	~PseudoRandom(){};  // Destructor
	
	/*****Sample from Various Distribution Functions *****/
	// This function returns a Gaussian random variate, with mean mu and standard deviation sigma
	double gaussian(const double mu, const double sigma);
	// This function returns a random number from a normal truncated to (left,right) interval with mean = mu, variance = sigma
	double normt_rnd(double mu,double sigma,double left,double right);
	// This function returns a random number from a left-truncated normal distribution, with mean = mu, variance = sigma
	double normlt_rnd(double mu,double sigma,double left);
	// This function returns a random number from a right-truncated normal distribution, with mean = mu, variance = sigma
	double normrt_rnd(double mu,double sigma,double right);
	// This function provides random variates from the upper tail of a Gaussian distribution with standard deviation sigma. The values returned are larger than the lower limit a, which must be positive
	double gaussian_tail(const double a, const double sigma);
	// This function generates a pair of correlated gaussian variates, with mean zero, correlation coefficient rho and standard deviations sigma_x and sigma_y in the x and y directions
	void bivariate_gaussian(double sigma_x, double sigma_y, double rho, double *x, double *y);
	// This function returns a random variate from the exponential distribution with mean mu
	double exponential(const double mu);
	// This function returns a random variate from the Laplace distribution with width a
	double laplace(const double a);
	// This function returns a random variate from the Cauchy distribution with scale parameter a
	double cauchy(const double a);
	// This function returns a random variate from the Rayleigh distribution with scale parameter sigma
	double rayleigh(const double sigma);
	// This function returns a random variate from the tail of the Rayleigh distribution with scale parameter sigma and a lower limit of a
	double rayleigh_tail(const double a, const double sigma);
	// This function returns a random variate from the gamma distribution
	double gamma(const double a, const double b);
	// This function returns a random variate from the uniform(flat) distribution from a to b
	double uniform(const double a, const double b);
	// This function returns a random variate from the lognormal distribution
	double lognormal(const double zeta, const double sigma);
	// This function returns a random variate from the chi-squared distribution with nu degrees of freedom
	double chisq(const double nu);
	// This function returns a random variate from the F-distribution with degrees of freedom nu1 and nu2
	double fdist(const double nu1, const double nu2);
	// This function returns a random variate from the t-distribution
	double tdist (const double nu);
	// This function returns a random variate from the beta distribution
	double beta(const double a, const double b);
	// This function returns a random variate from the logistic distribution
	double logistic(const double a);
	// This function returns a random variate from the Pareto distribution of order a
	double pareto(double a, const double b);
	// This function returns a random variate from the Weibull distribution
	double weibull(const double a, const double b);
	// This function returns a random variate from the Type-1 Gumbel distribution
	double gumbel1(const double a, const double b);
	// This function returns a random variate from the Type-2 Gumbel distribution
	double gumbel2(const double a, const double b);
	// This function returns an array of K random variates from a Dirichlet distribution of order K-1
	void dirichlet(const unsigned int K, const double alpha[], double theta[]);
	// This function returns a random integer from the Poisson distribution with mean mu
	unsigned int poisson(double mu);
	// This function returns either 0 or 1, the result of a Bernoulli trial with probability p
	unsigned int bernoulli(double p);
	// This function returns a random integer from the binomial distribution, the number of successes in n independent trials with probability p
	unsigned int binomial(double p, unsigned int n);
	// This function returns an array of K random variates from a multinomial distribution
	void multinomial(const unsigned int K, const unsigned int n, const double p[], unsigned int d[]);
	// This function returns a random integer from the negative binomial distribution, the number of failures occurring before n successes in independent trials with probability p of success
	unsigned int negative_binomial(double p, double n);
	// This function returns a random integer from the Pascal distribution. The Pascal distribution is simply a negative binomial distribution with an integer value of n
	unsigned int pascal_(double p, unsigned int n);
	// This function returns a random integer from the geometric distribution, the number of independent trials with probability p until the first success
	unsigned int geometric(const double p);
	// This function returns a random integer from the hypergeometric distribution
	unsigned int hypergeometric(unsigned int n1, unsigned int n2, unsigned int t);
	// This function returns a random integer from the logarithmic distribution
	unsigned int logarithmic(const double p);
	// This function returns a matrix associated with multivariate normal distribution with mean mu and covariance matrix sigma, cases is the number of samples
	Matrix<double> mnormrnd(int cases, Matrix<double> & mu, Matrix<double> & sigma);
	// This function returns a matrix associated with truncated multivariate normal distribution with mean mu and covariance matrix sigma, n is the number of samples
	Matrix<double> tnorm_rnd(int n, Matrix<double> & amu, Matrix<double> & sigma, Matrix<double> & a, 
	                         Matrix<double> & b,Matrix<double> & la, Matrix<double> & lb, Matrix<double> & d,Matrix<int> & kstep);
	
	/*****Calculate the Probability Density Functions(PDF) *****/
	// This function computes the probability density p(x) at x for a Gaussian distribution with mean mu and standard deviation sigma
	double gaussian_pdf(const double x, const double mu, const double sigma);
	// This function computes the probability density p(x,y) at (x,y) for a bivariate gaussian distribution with standard deviations sigma_x, sigma_y and correlation coefficient rho, using the formula given above
	double bivariate_gaussian_pdf(const double x, const double y, const double sigma_x, const double sigma_y, const double rho);
	// This function computes the probability density p(x) at x for an exponential distribution with mean mu
	double exponential_pdf(const double x, const double mu);
	// This function computes the probability density p(x) at x for a uniform distribution from a to b
	double flat_pdf(double x, const double a, const double b);
	/*****Calculate the Cumulative Distribution Functions(PDF) *****/
	

private:
	/***** Subfunctions used by several functions *****/
	double gamma_int(const unsigned int a);
	double gamma_large(const double a);
	double gamma_frac(const double a);
	double normal_01_cdf ( double x );
	double normal_01_cdf_inv ( double p );
	double d_huge ( void );
	double dpoly_value ( int n, double a[], double x );	
};
		
#endif /*PseudoRandom_H_*/