<div><br></div>You'll need some patience to get non-zeros, especially for k=1e-5<div><br></div><div><div>In [84]: np.sum(np.random.gamma(1e-5,size=1000000)!=0.0)</div><div>Out[84]: 7259</div><div>that's less than 1%. For k=1e-4 it's ~7%</div>
<div><br></div><div>Val</div><div><br></div><br><div class="gmail_quote">On Mon, May 28, 2012 at 10:33 PM, Uri Laserson <span dir="ltr"><<a href="mailto:uri.laserson@gmail.com" target="_blank">uri.laserson@gmail.com</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">I am trying to sample from a Dirichlet distribution, where some of the shape parameters are very small.  To do so, the algorithm samples each component individually from a Gamma(k,1) distribution where k is the shape parameter for that component of the Dirichlet.  In principle, this should always return a positive number (as the Dirichlet is defined).  However, if k is very small, it will return zero:<div>


<br></div><div><div>In [157]: np.random.gamma(1e-1)</div><div>Out[157]: 4.863866491339177e-06</div><div><br></div><div>In [158]: np.random.gamma(1e-2)</div><div>Out[158]: 2.424451829710714e-57</div><div><br></div><div>In [159]: np.random.gamma(1e-3)</div>


<div>Out[159]: 5.1909861689757784e-197</div><div><br></div><div>In [160]: np.random.gamma(1e-4)</div><div>Out[160]: 0.0</div><div><br></div><div>In [161]: np.random.gamma(1e-5)</div><div>Out[161]: 0.0</div><div><br></div>


<div>What is the best way to deal with this?</div><div><br></div><div>Thanks!</div><div>Uri</div><div><br></div><br clear="all"><div><br>...................................................................................</div>
<span class="HOEnZb"><font color="#888888">

Uri Laserson<br>Graduate Student, Biomedical Engineering<br>Harvard-MIT Division of Health Sciences and Technology<br>M <a href="tel:%2B1%20917%20742%208019" value="+19177428019" target="_blank">+1 917 742 8019</a><br><a href="mailto:laserson@mit.edu" target="_blank">laserson@mit.edu</a><br>

</font></span></div>
<br>_______________________________________________<br>
NumPy-Discussion mailing list<br>
<a href="mailto:NumPy-Discussion@scipy.org">NumPy-Discussion@scipy.org</a><br>
<a href="http://mail.scipy.org/mailman/listinfo/numpy-discussion" target="_blank">http://mail.scipy.org/mailman/listinfo/numpy-discussion</a><br>
<br></blockquote></div><br></div>