2

I am working on implementation of a machine learning method that in part of the algorithm I need to calculate the value of $\Gamma (\alpha) / \Gamma (\beta) $. $\alpha$ and $\beta$ are quite large numbers (i.e. bigger than 200) and it causes the python $gamma$ function to overflow. However, as the difference of $\alpha$ and $\beta$ is relatively small (e.g. $|\alpha-\beta|<5$), the final result is not such a big number and can be used for later purposes. So, I am trying to calculate (or approximate) the value of $\Gamma (\alpha) / \Gamma (\beta) $ without going through the calculation of $\Gamma (\alpha)$ and $\Gamma (\beta)$ directly. If $\alpha$ and $\beta$ were integers, the result would be simple equal to $\alpha . \alpha+2. \alpha+3... \beta-1$, But I can not imagine how this formula will be changed if we let $\alpha$ and $\beta$ to be real numbers.

2 Answers2

4

I think that a good solution would be Stirling approximation that is to say $$\log(\Gamma(x))=x (\log (x)-1)+\frac{1}{2} \left(-\log \left({x}\right)+\log (2 \pi )\right)+\frac{1}{12 x}+O\left(\frac{1}{x^3}\right)$$ Now, consider $$y=\frac{\Gamma(\alpha)}{\Gamma(\beta)}\implies \log(y)=\log(\Gamma(\alpha))-\log(\Gamma(\beta))$$ Apply the formula (even with more terms) and use later $y=e^{\log(y)}$.

You are then able to control overflows and underflows if required.

By the way, why not to use in Python function lgamma(x) ?

  • Thank you. I was not aware of lgamma function. It seems that it solves the problem in addition to using the Stirling approximation directly. – CoderInNetwork Dec 25 '16 at 09:48
  • 1
    @CoderInNetwork. Since I do not use Python, I did not know either. I just googled for mathematical functions Python. Using lgamma(x) will be much better than Stirling (I am sure). Cheers. – Claude Leibovici Dec 25 '16 at 15:10
3

From this answer we have Gautschi's inequality:

$$x^{1-s}<\frac{\Gamma(x+1)}{\Gamma(x+s)}<(x+1)^{1-s},\quad x>0,0<s<1$$

We can combine this with the functional equation for the gamma function: $$\Gamma(z+1) = z\Gamma(z)$$ which holds for all $z$ where $\Gamma$ is defined.

Now, assume that $\alpha\leq\beta$, and that $|\beta-\alpha|\leq 5$. Now, let $k$ be the greatest integer such that $\alpha+k\leq\beta$ (so clearly, $k\in\{0,1,2,3,4,5\}$. It follows that $\alpha+k+s = \beta$ where $0\leq s<1$.

Now, we have that $\Gamma(\beta) = \Gamma(\alpha+k+s)$, and we can apply the functional equation $k+1$ times to get that $\Gamma(\beta) =\prod_{i = 0}^{k+1} (\alpha+k+s)\Gamma(\alpha-1+s)$. It follows that: $$\frac{\Gamma(\alpha)}{\Gamma(\beta)} = \frac{\Gamma(\alpha)}{\Gamma(\alpha+k+s)} = \frac{\Gamma(\alpha)}{\prod_{i = 0}^k (\alpha+k+s)\Gamma(\alpha-1+s)}$$ Now, if we substitue in $\alpha-1 = x$, we see that: $$\frac{\Gamma(\alpha)}{\Gamma(\beta)} = \frac{1}{\prod_{i = 0}^{k+1}((\alpha+k+s)}\frac{\Gamma(x+1)}{\Gamma(x+s)}$$ This part in the front is just some constant $c$, so we can rewrite this as: $$\frac{\Gamma(\alpha)}{\Gamma(\beta)} = c\frac{\Gamma((\alpha-1)+1)}{\Gamma((\alpha-1)+s)}$$ From Gautschi's inequality, we have that: $$(\alpha-1)^{1-s}<\frac{\Gamma(\alpha-1+1)}{\Gamma(\alpha-1+s)}<\alpha^{1-s}$$ Now, as you said that $\alpha\gg 5$, we have that $c\geq 0$, so we can multiply through by $c$ to get that: $$c(\alpha-1)^{1-s}<\frac{\Gamma(\alpha)}{\Gamma(\beta)}<c\alpha^{1-s}$$ The above calculation implicitly assume that $s\neq 0,1$, but if this is true then we can get the result just from the functional equation for the gamma. Moreover, $0<1-s<1$, so while $\alpha$ may be large, $\alpha^{1-s}$ is hopefully manageable.

This is all to try to give you another option.

  • Thanks for your answer. up-voted. – CoderInNetwork Dec 25 '16 at 09:46
  • You've written “$\Gamma(\beta) =\prod_{i = 0}^{k+1} (\alpha+k+s)\Gamma(\alpha-1+s)$”, where $\beta = \alpha + k + s$. But, why does $i$ not show up in the product? My analysis suggests that $\Gamma(\beta) = \Gamma(\alpha - 1 + s) \prod_{i=0}^k(\alpha - 1 + s + i)$ and I don't think that quite agrees? – Ahmed Fasih Nov 19 '18 at 18:56