In [1]:

    

import numpy as np


    

In [131]:

    

bin_count = 50
bins = np.array([0]*bin_count)

bins[5] += 10
bins[6] += 1


bins[46] += 1

bin_values = np.array([i for i in range(0,bin_count)])


    

In [132]:

    

possibility_distribution = [(bin_values*np.random.dirichlet(bins, bins.sum())/bins.sum()).sum() for i in range(5000)]


    

In [133]:

    

%matplotlib inline
import matplotlib.pyplot as plt

plt.xlim(0,bin_count)
plt.hist(possibility_distribution, bins=50)
''


    

    
Out[133]:





''






    

In [13]:

    

(bin_values*np.random.dirichlet(bins, 2)).mean()


    

    
Out[13]:





0.11303761322863445

In [54]:

    

bins = np.array([1]*20)
bins[10] += 1
bins[15] += 1
observations = 1

np.random.dirichlet(bins, 1)


    

    
Out[54]:





array([[ 0.17338139,  0.03455451,  0.02325679,  0.02124659,  0.02573916,
         0.14553389,  0.00479338,  0.04005475,  0.04634759,  0.01063245,
         0.10830318,  0.03707728,  0.03097544,  0.0554245 ,  0.00107638,
         0.11514135,  0.01751338,  0.05830946,  0.02503389,  0.02560464]])

In [30]:

    

(bin_values*np.random.dirichlet(bins, 2)/2.0).sum()


    

    
Out[30]:





12.944118435272863

In [20]:

    

np.random.multinomial?


    

In [53]:

    

np.random.dirichlet?


    

In [136]:

    

def custom_dirichlet(observations):
  total_observations = len(observations) + observations.sum()
  sample = np.array([np.random.beta(1+total_observations-obs, 1+obs) for obs in observations])
  return sample


    

In [137]:

    

custom_dirichlet(bins)


    

    




---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-137-922c2382d813> in <module>()
----> 1 custom_dirichlet(bins)

<ipython-input-136-f6cb19ef9a1f> in custom_dirichlet(observations)
      1 def custom_dirichlet(observations):
      2   total_observations = len(observations) + observations.sum()
----> 3   sample = np.array([np.random.beta(1+observations-obs, 1+obs) for obs in observations])
      4   return sample

<ipython-input-136-f6cb19ef9a1f> in <listcomp>(.0)
      1 def custom_dirichlet(observations):
      2   total_observations = len(observations) + observations.sum()
----> 3   sample = np.array([np.random.beta(1+observations-obs, 1+obs) for obs in observations])
      4   return sample

mtrand.pyx in mtrand.RandomState.beta (numpy/random/mtrand/mtrand.c:12967)()

ValueError: a <= 0

In [ ]: