Hi,

I encountered a bug when using np.random.dirichlet with small alpha parameters. Call and traceback are below.

ZeroDivisionError                         Traceback (most recent call last)
<ipython-input-86-73c2067e20c1> in <module>()
----> 1 np.random.dirichlet([0.0001, 0.0, 0.0001])

mtrand.pyx in mtrand.RandomState.dirichlet (numpy/random/mtrand/mtrand.c:24477)()

mtrand.pyx in mtrand.RandomState.dirichlet (numpy/random/mtrand/mtrand.c:24387)()

ZeroDivisionError: float division

I am using numpy-1.9.1.

I believe this is a floating point issue, the distribution has almost all of its mass very close to either (1, 0, 0) or (0, 0, 1) .

The 'float division' error already occurs for larger values<1, e.g. 0.001.

It is likely that this occurs because of the Dirichlet distribution is usually sampled via the Gamma distribution followed by normalization. If all values returned from Gamma sampling are zero than a float division error occurs.

In addition

np.random.beta(0.0001, 0.0001)

produces 'nan' most of the time, while it should be alternating 'almost always' between (1, 0) and
(0, 1)

python scipy.special.betainc(0.0001, 0.0001, 1e-50) = 0.494..

It might not be able to fix that in the current algorithmic framework but maybe it is possible to discourage/prevent users from supplying too small parameters.

Wow, this wasn't supposed to become such a long post. Thanks to anyone reading/considering this issue.