PR: Add logaddexp to the spec

[steff456]

This PR

• adds specification for logaddexp, which calculates the logarithm of the sum of exponentiations of the inputs.
• follows the conventions set forth by existing element-wise functions.

[rgommers]

Thanks @steff456

All the libraries use the alias logaddexp. Went with logsumexp is more clearer and easy to read/say.

The specification looks good, but I'm not sure about this name change. If all libraries name it logaddexp, and that's also the more logical name (because it's not a reduction, it's an elementwise operation), why would we want to rename it?

[asmeurer]

MATLAB uses logsumexp, but it operates on a vector (see https://nhigham.com/2021/01/05/what-is-the-log-sum-exp-function/), and it's also the "mathematical" name (https://en.wikipedia.org/wiki/LogSumExp). Presumably the fact that it operates on two arguments is why it was changed to logaddexp.

Is reduction part of the spec? I am surprised that logaddexp is defined as a two argument function rather than a function on a vector, but I suppose that's because you can just use logaddexp.reduce since it's associative (and as far as I can tell, computing it like this is numerically stable).

[rgommers]

Is reduction part of the spec?

Reductions as in sum, std, etc. yes. A .reduce method, definitely not.

I am surprised that logaddexp is defined as a two argument function rather than a function on a vector, but I suppose that's because you can just use logaddexp.reduce since it's associative (and as far as I can tell, computing it like this is numerically stable).

You're talking about a different function altogether. From your first link: "The log-sum-exp function takes as input a real n-vector x and returns the scalar"

From the np.logaddexp docs:

Logarithm of the sum of exponentiations of the inputs.

Calculates ``log(exp(x1) + exp(x2))``. This function is useful in
statistics where the calculated probabilities of events may be so small
as to exceed the range of normal floating point numbers.  In such cases
the logarithm of the calculated probability is stored. This function
allows adding probabilities stored in such a fashion.

It's not a reduction, as sum and your other logsumexp are:

>>> x = np.arange(5.)
>>> y = np.ones(5)
>>> np.logaddexp(x, y)
array([1.31326169, 1.69314718, 2.31326169, 3.12692801, 4.04858735])

[rgommers]

scipy.special on the other hand does have logsumexp:

>>> from scipy.special import logsumexp
>>> x = np.arange(5.)
>>> logsumexp(x)
4.451914395937593

[asmeurer]

I think you missed the point of what I said. logaddexp is associative. Its reduction (via reduce) is logsumexp.

>>> np.logaddexp.reduce(np.arange(5, dtype=np.float64))
4.451914395937593

Most references talk about the n-argument version (logsumexp), so I'm surprised that NumPy only has a two argument version, but I imagine that the reason is that you can easily get the n-argument version via reduce (but that isn't part of the spec).

[steff456]

Hi guys,

You are right! I didn't realize that logsumexp and logaddexp are different functions in Torch and TF. I'll revert the change of the name and I'll double check the signature.

[kgryte]

LGTM. Will add two minor edits, but otherwise looks good! Thanks @steff456!

[kgryte]

As this PR has been open a week, the proposal straightforward, and can be readily implemented based on existing standardized APIs, will merge.