Gamma function returns unstable value?

Question

Gamma function should not take any negative value as an argument. Look at the code below where strange thing happens. Is this some problem with R?

I was using function optim to optimize some function containing:

gamma(sum(alpha))

with respect to alpha. R returns negative alpha.

> gamma(sum(alpha))
[1] 3.753+14
>sum(alpha)
[1] -3
gamma(-3)
[1] NaN
Warning message:
In gamma(-3) NaN's produced.

Can somebody explain? Or any suggestion for the optimization?

Thanks!

-3 is a pole, and thus gamma function as x tends to -3 will be infinity. for you to note that `sum(alpha)` is not `-3` just do simply `sum(alpha)==-3` it will give you false. Since `sum(alpha)` is technically close to -3 but is not -3, `identical(sum(alpha),-3)` — Onyambu, Feb 12 '18 at 06:43
It appears that everyone now needs to read the `?gamma` help page. — IRTFM, Feb 12 '18 at 06:46
From the help: `The gamma function is defined by (Abramowitz and Stegun section 6.1.1, page 255) Γ(x) = integral_0^Inf t^(x-1) exp(-t) dt for all real x except zero and negative integers (when NaN is returned). ` — kangaroo_cliff, Feb 12 '18 at 06:48
In R -3 is not of integer type. Looks like `sum(alpha)` is near `-3L-0.0000000000000004` — IRTFM, Feb 12 '18 at 06:51
Thanks! So Gamma function is defined for negative numbers expect negative integers. — Xinming, Feb 12 '18 at 07:23
Yes. We can show that `Gamma(x) = 1/x Gamma(x)`. If `x = -3/2`, then `Gamma(-3/2) = -2/3 * Gamma(-1/2) = -2/3 * -2/1 * Gamma(1/2) = 4/3 * sqrt(pi)`. It fails for negative integers because `Gamma(0)` isn't defined. (grad school notes: http://www.suchanutter.net/ItCanBeShown/gamma-function.html) — Benjamin, Feb 12 '18 at 12:08

Artem · Accepted Answer · 2018-10-05T18:03:14.187

Gamma function is "not defined" at negative integer argument values so R returns Not a Number (NaN). The reason of the "strange" behaviour is decimal representation of numbers in R. In case the number differs from the nearest integer not very much, R rounds it during printing (in fact when you type alpha, R is calling for print(alpha). Please see the examples of such a behaviour below.

gamma(-3)
# [1] NaN
# Warning message:
#  In gamma(-3) : NaNs produced
x <- -c(1, 2, 3) / 2 - 1e-15
x
# [1] -0.5 -1.0 -1.5
sum(x)
# [1] -3
gamma(sum(x))
# [1] 5.361428e+13
curve(gamma, xlim = c(-3.5, -2.5))

Please see a graph below which explains the behaviour of gamma-function near negative integers:

Gamma function returns unstable value?

1 Answers1