Geek Answers Handbook

Beta Random Random Number Generator

I need c or C ++ source code for a function like betarand(a,b) that creates a random number with a beta distribution. I know that I can use the boost library, but I'm going to port it for the CUDA architecture, so I need code. Can someone help me?
Meanwhile, I have betapdf (beta probability density function). But I do not know how to use it to create random numbers :).

+6
source share
4 answers

The C ++ 11 random number library does not provide a beta distribution. However, the beta distribution can be modeled in terms of the two gamma distributions that the library provides. I have executed beta_distribution in terms of std::gamma_distribution for you. As far as I can tell, it fully complies with the requirements for the distribution of random numbers.

 #include <iostream> #include <sstream> #include <string> #include <random> namespace sftrabbit { template <typename RealType = double> class beta_distribution { public: typedef RealType result_type; class param_type { public: typedef beta_distribution distribution_type; explicit param_type(RealType a = 2.0, RealType b = 2.0) : a_param(a), b_param(b) { } RealType a() const { return a_param; } RealType b() const { return b_param; } bool operator==(const param_type& other) const { return (a_param == other.a_param && b_param == other.b_param); } bool operator!=(const param_type& other) const { return !(*this == other); } private: RealType a_param, b_param; }; explicit beta_distribution(RealType a = 2.0, RealType b = 2.0) : a_gamma(a), b_gamma(b) { } explicit beta_distribution(const param_type& param) : a_gamma(param.a()), b_gamma(param.b()) { } void reset() { } param_type param() const { return param_type(a(), b()); } void param(const param_type& param) { a_gamma = gamma_dist_type(param.a()); b_gamma = gamma_dist_type(param.b()); } template <typename URNG> result_type operator()(URNG& engine) { return generate(engine, a_gamma, b_gamma); } template <typename URNG> result_type operator()(URNG& engine, const param_type& param) { gamma_dist_type a_param_gamma(param.a()), b_param_gamma(param.b()); return generate(engine, a_param_gamma, b_param_gamma); } result_type min() const { return 0.0; } result_type max() const { return 1.0; } result_type a() const { return a_gamma.alpha(); } result_type b() const { return b_gamma.alpha(); } bool operator==(const beta_distribution<result_type>& other) const { return (param() == other.param() && a_gamma == other.a_gamma && b_gamma == other.b_gamma); } bool operator!=(const beta_distribution<result_type>& other) const { return !(*this == other); } private: typedef std::gamma_distribution<result_type> gamma_dist_type; gamma_dist_type a_gamma, b_gamma; template <typename URNG> result_type generate(URNG& engine, gamma_dist_type& x_gamma, gamma_dist_type& y_gamma) { result_type x = x_gamma(engine); return x / (x + y_gamma(engine)); } }; template <typename CharT, typename RealType> std::basic_ostream<CharT>& operator<<(std::basic_ostream<CharT>& os, const beta_distribution<RealType>& beta) { os << "~Beta(" << beta.a() << "," << beta.b() << ")"; return os; } template <typename CharT, typename RealType> std::basic_istream<CharT>& operator>>(std::basic_istream<CharT>& is, beta_distribution<RealType>& beta) { std::string str; RealType a, b; if (std::getline(is, str, '(') && str == "~Beta" && is >> a && is.get() == ',' && is >> b && is.get() == ')') { beta = beta_distribution<RealType>(a, b); } else { is.setstate(std::ios::failbit); } return is; } }

Use it like this:

 std::random_device rd; std::mt19937 gen(rd()); sftrabbit::beta_distribution<> beta(2, 2); for (int i = 0; i < 10000; i++) { std::cout << beta(gen) << std::endl; }
+14
source

Perhaps you can use the code that gsl uses to create random numbers with beta distribution. They use a slightly strange way to create them, since you need to pass a random number generator to a function, but, of course, you can get what you need.

Here is the documentation and the web page

+2
source

Boost Reverse Incomplete Beta is another quick (and easy) way to simulate a beta.

 #include <random> #include <boost/math/special_functions/beta.hpp> template<typename URNG> double beta_sample(URNG& engine, double a, double b) { static std::uniform_real_distribution<double> unif(0,1); double p = unif(engine); return boost::math::ibeta_inv(a, b, p); // Use Boost policies if it not fast enough }
+1
source

Check out NumPy random number generator implementation: NumPy distribution source

They are implemented in C and are very fast.

0
source

More articles:


All Articles