2

I cannot find much documentation of the Boost version of discrete_distribution. After much Google searching, I still can't even find a list of methods that this class has, and whether any of them function to re-assign the probabilities.

In my case, I am writing an evolutionary dynamics algorithm. At each time step, members of the population can be randomly selected to either die or reproduce. Because of this, the total number of entries within my discrete distribution will change almost every iteration.

I want to have a single object that I define before the simulation starts, called gillespie_dist (the discrete distribution governing this Gillespie algorithm). But I want, at the end of each iteration, to potentially change specific values and/or add new values to gillespie_dist and specifically do not want to create new instances of the discrete_distribution every iteration.

What is a good approach for this. Are there methods for pushing a new value onto the discrete_distribution object, methods for changing a distribution's value at a particular index, or better yet, somehow "re-initializing" the entire distribution using the vector iterator idea mentioned here?

Community
  • 1
  • 1
ely
  • 74,674
  • 34
  • 147
  • 228
  • `std::discrete_distribution` can be entirely re-initialized with the [`param`](http://en.cppreference.com/w/cpp/numeric/random/discrete_distribution/param) member function, but that would hardly be less expensive than creating a new `discrete_distribution`. I don't think these classes are meant to support small updates efficiently. – leftaroundabout Mar 31 '12 at 01:26
  • As is implied by previous comment, std::discrete_distribution is a standard feature of C++ as of last year: http://en.cppreference.com/w/cpp/numeric/random/discrete_distribution – Andrew Tomazos Mar 31 '12 at 01:49
  • Noted, but as in my comment to the answer below, I think this has been a version issue. I get `error: ‘discrete_distribution’ is not a member of ‘std’` when I try to use `std::discrete_distribution`. I'll be updating shortly and things should be smoother then. – ely Mar 31 '12 at 01:52
  • @EMS: What compiler and version are are you using? – Andrew Tomazos Mar 31 '12 at 01:59
  • g++ (Ubuntu/Linaro 4.4.5-15ubuntu1) 4.4.5 – ely Mar 31 '12 at 02:00
  • Not sure if it is in that release, but to use it you need to turn C++11 features on with "-std=gnu++0x" compiler switch. But have another look at my answer below, you might see how easy it is to implement it yourself. Also if you upgrade to Oneiric (Ubuntu 11.10) or Precise beta (Ubuntu 12.04 beta), than you get gcc 4.6 which has much better C++11 support. – Andrew Tomazos Mar 31 '12 at 02:14
  • As I mentioned in a reply to an answer below, the `param` function does not appear to work to reassign things with a std::vector or any similar resizable array. It also doesn't accept atd::vector iterators to do the re-initializing. The whole point of this idea is to have a dynamically resizable discrete distribution, so it looks like the only thing to do is write my own sampling code that works with a std::vector that contains probabilities. – ely Apr 01 '12 at 20:19

2 Answers2

1

I looked into the code of the gcc libstdc++ 4.7 implementation of std::discrete_distribution.

The weights are stored as a vector<double> in a private member. There is no access to its resize method in the public interface.

I'll try and dig out the implementation of its operator() (which is in the cpp it looks like), it should be no trouble to roll your own.

Here is the main action, and my explanation following:

  template<typename _IntType>
    void
    discrete_distribution<_IntType>::param_type::
    _M_initialize()
    {
      if (_M_prob.size() < 2)
        {
          _M_prob.clear();
          return;
        }

      const double __sum = std::accumulate(_M_prob.begin(),
                                           _M_prob.end(), 0.0);
      // Now normalize the probabilites.
      __detail::__transform(_M_prob.begin(), _M_prob.end(), _M_prob.begin(),
                          std::bind2nd(std::divides<double>(), __sum));
      // Accumulate partial sums.
      _M_cp.reserve(_M_prob.size());
      std::partial_sum(_M_prob.begin(), _M_prob.end(),
                       std::back_inserter(_M_cp));
      // Make sure the last cumulative probability is one.
      _M_cp[_M_cp.size() - 1] = 1.0;
    }

  template<typename _IntType>
    template<typename _UniformRandomNumberGenerator>
      typename discrete_distribution<_IntType>::result_type
      discrete_distribution<_IntType>::
      operator()(_UniformRandomNumberGenerator& __urng,
                 const param_type& __param)
      {
        if (__param._M_cp.empty())
          return result_type(0);

        __detail::_Adaptor<_UniformRandomNumberGenerator, double>
          __aurng(__urng);

        const double __p = __aurng();
        auto __pos = std::lower_bound(__param._M_cp.begin(),
                                      __param._M_cp.end(), __p);

        return __pos - __param._M_cp.begin();
      }

So basically it calculates an auxilary vector _M_cp at initialization time which is essentially a discete cummulative density function of the weights. So producing a sample just means generating a uniform random variable and searching for the first occurence of that in the cummulative distribution (this is the lower_bound call above), returning its index.

eg, if the weights vector is:

{ 1, 2, 1, 3 }

then the cp is calculated as:

{ 1, 1+2, 1+2+1, 1+2+1+3 }
=
{ 1, 3, 4, 7 }

so I uniformly select from 0..6 and get 4, so I pick the third one.

Andrew Tomazos
  • 66,139
  • 40
  • 186
  • 319
0

After much Google searching, I still can't even find a list of methods that this class has, and whether any of them function to re-assign the probabilities.

http://www.boost.org/doc/html/boost/random/discrete_distribution.html

and

void param(const param_type & param);

Sets the parameters of the distribution.

ildjarn
  • 62,044
  • 9
  • 127
  • 211
DRVic
  • 2,481
  • 1
  • 15
  • 22
  • I see. The problem was that my system has Boost 1.42, and it appears that this stuff doesn't exist in Boost::random in that version. I will upgrade to the current Boost version and then the documentation you pointed me to will apply and should solve the problem. Once I get that working, I'll be sure to accept the answer. Thanks. – ely Mar 31 '12 at 01:37
  • I cannot get this method to work. Is there documentation on what is accepted as valid for 'param' type? I have tried `std::vector` which is what seems right from examples, but it gives errors that it cannot find a `param()` method with that type as the argument. – ely Mar 31 '12 at 04:06
  • 1
    It actually wants a param_type which can be initialized in various ways, but vector is not one of them. An array initialized as { 1.0, 2.0, 1.3 } would work, but that's not very flexible. There's also the iterator, iterator form which might work for you. – DRVic Mar 31 '12 at 12:07
  • I get the same error when I try to pass it two std::vector iterators. It says it can't find a method with those types as inputs. The issue is that I need this for dynamic resizing, and so having to form a bounded-memory array ruins all the gain I would get from using this method anyway. I need some discrete_distribution object which can accept dynamically resizable containers and update its probabilities to match them, whether this means adding to or reducing the previous size of the distribution. This doesn't appear to do that, unless I just make some massive regular array at the start – ely Mar 31 '12 at 17:29