Why is numpy.mean infinite while values are finite?

Question

I was playing around with Hypothesis when I came across the following

>>> x = np.array([8.988465674311579e+307, 8.98846567431158e+307])
>>> x
array([  8.98846567e+307,   8.98846567e+307])
>>> x.mean()
inf

Does this happen because the sum of the elements is inf?

>>> x.sum()
inf

If so, why does numpy not divide first, sum after?

>>> (x/len(x)).sum()
8.9884656743115795e+307

Answer 1

x.mean is a builtin, so I'd have to look at the github to see exactly how it's coded. But the straight forward sum followed by divide is most likely. But it does handle variations like axis and dtype.

My guess is that your case does not come up often enough to merit special treatment. Routinely dividing by size is not a good idea, since that can lead to a loss of precision. And it involves an extra iteration through the array. The alternative is to do the sum(s), and then check if any are inf or overflow values, and repeat the task with the scaled values. That too incurs an extra iteration.

But you as the user always have the option of scaling your values appropriately. Coding your own mean with the scaling followed by sum is easy. If you hit this inf problem with mean you probably will hit with other calculations with those numbers - anything that involves summing - weighted average, dot, std, etc. So you might as well scale the values once, rather than expect each calculation to do it for you.

Why does "numpy.mean" return 'inf'? suggests that this is more likely a problem when you are deliberately limiting the dtype (e.g. with dtype=np.float16).

In [361]: x=np.array([6550]*1000,dtype=np.float16)
In [362]: x.mean()
Out[362]: inf
In [363]: x.mean(dtype=float)
Out[363]: 6552.0

np.mean docs talk about the choice of dtype and potential loss of precision. They could have just as well mentioned the overflow issue.