# Caculating the mean of vector accurately

I am having trouble with calculating a mean of vector with sufficient accuracy. My current solution which works but it quite slow and has unpredictable performance:

mean_sum = mean = math.fsum(values) / n
values = values - mean
mean = math.fsum(values) / n
while np.abs(mean) > np.finfo(np.float32).eps:
values = values - mean
mean_sum += mean
mean = math.fsum(values) / n
return mean_sum, values


Is there a better way?

EDIT: All the values are non-negative. The accuracy required is defined by

np.finfo(np.float32).eps


as the mean of the resulting series. I have edited the code a bit to make it more clear.

• Is there a reason that you seem to be calculating a running average, i.e., you divide by n inside your loop, rather than just the mean? Look at the "pseudocode" link in my answer and not how only the sums are accumulated with the loop. If you really do need to access the current mean ass you sum your vector, do that in addition to a variable for accumulating. – horchler Nov 3 '13 at 20:10
• n is simply the length of the vector. I need to have the centered vector and the its original mean. If I only calculate it once, the resulting vector has a mean with a magnitude of around -4 which is unacceptable. Perhaps the problem is not with the mean but rather with the division or subtraction? The vector length is about 8000. I checked the Kahan Summation from your answer, it looks interesting but I think fsum is already doing a good job and something similar as @k20 said. Are you sure the problem is with the summation? – Sonia Nov 3 '13 at 20:20

How accurately? Is there a reason you're not just using numpy.mean? Is that not sufficient? If you need to compensate for floating-point error, you can try using Kahan summation. Here's some pseudocode (PDF).

For a review of floating-point summation techniques see: N.J. Higham, SIAM, 1993.

• The OP is already using math.fsum() which supposedly guarantees rounding to within floating point precision svn.python.org/projects/python/trunk/Modules/mathmodule.c – k20 Nov 3 '13 at 19:58
• @k20: I noticed that. I think the true issue may be the use of a while loop in conjunction with calculating a running average rather than a sum that is divided once at the end. – horchler Nov 3 '13 at 20:12

If you're concerned about floating point round-off, then have a loop over all elements, add the positive elements to one variable and the negative elements to another, and then add the two together to get the sum over all elements. Divide the result by the number of elements and you have the best mean value you can likely get.

• Can't you do better by sorting each of the two lists from smallest to largest before adding down each list? – Bill Barth Nov 3 '13 at 20:02
• Yes. And if you're concerned with the $O(N \log N)$ cost, it is enough to sort into buckets of elements of roughly the same size. This can be done in $O(N)$ effort. – Wolfgang Bangerth Nov 3 '13 at 22:36

I see that you are using Python, so you could use a high precision library like mpmath.