9 ms·
How to average floating point numbers (2015)
- deleted 5y ago[deleted]
- eesmith 5y agoSwitching to Python instead of Perl, use the 'statistics' module, https://docs.python.org/3/library/statistics.html https://docs.python.org/3/library/statistics.html : >>> data = (1_000_000_000.1, 1.1) * 50_000 >>> sum(data) / len(data) 500000000.60091573 >>> sum(sorted(data)) / len(data) 500000000.6004578 >>> sum(sorted(data, reverse=True)) / len(data) 500000000.6012391 >>> import statistics >>> statistics.mean(data) 500000000.6 and >>> X = [10_000_000_000 + k for k in (4, 7, 13, 16)] >>> statistics.variance(X) 30
- pdkl95 5y agoRuby seems to handle this out of the box: >> RUBY_VERSION ⇒ "2.7.2" >> data = [1_000_000_000.1, 1.1] * 50_000 ; nil ⇒ nil >> data.sum ⇒ 50000000060000.0 >> data.sum / data.length ⇒ 500000000.6 >> data.sort.sum / data.length ⇒ 500000000.6 >> data.sort.reverse.sum / data.length ⇒ 500000000.6 I know that ruby will auto-promote an Integer into a Bignum. It seems like ruby will also auto-promote other types? I should research this...
- dmurray 5y agoFWIW, Ruby the language doesn't have Bignum any more, though, just Integer which handles arbitrarily large numbers. Under the hood, there are some optimizations for small integers. https://bugs.ruby-lang.org/issues/12005 https://bugs.ruby-lang.org/issues/12005
- chrisseaton 5y agoBignum is still there, working in almost exactly the same way, it's not exposed to the user anymore.
- aardvark179 5y agoNope. The sum method is implemented in C, and as soon as it sees a float value in the array it switches to Kahan-Babuska summation. We don’t do that in TruffleRuby and so don’t produce quite such an accurate result.
- ska 5y agoThis is a nice thing about Common Lisp, not just integer promotion to bignums, but ratios too (works in complex, also).
- trinovantes 5y agoI guess the question is does it still use floats internally or use something like BigNum and handle the precision issue (and thus more memory/cpu)
- _delirium 5y agoIt looks like the second. For many of its operations, the statistics module first converts floats to integer ratios (integers are bignums in Python), does intermediate operations using integer arithmetic, and converts back to floats at the end. Here's the internal _sum function for example, which you can see maps _exact_ratio over the data before summing: https://github.com/python/cpython/blob/3.9/Lib/statistics.py#L123 https://github.com/python/cpython/blob/3.9/Lib/statistics.py...
- eesmith 5y agoThere's also fmean() which is for floats: >>> data = (1_000_000_000.1, 1.1) * 50_000 >>> statistics.fmean(data) 500000000.6 Implementation at https://github.com/python/cpython/blob/6df926f1c46eb6db7b5dcd0227c6b532c78525c9/Modules/mathmodule.c#L1314 https://github.com/python/cpython/blob/6df926f1c46eb6db7b5dc... .
- eesmith 5y agoI just found that boltons supports a mean, in boltons.statsutils.mean . It appears to use the naive summation: >>> import boltons.statsutils >>> data = (1_000_000_000.1, 1.1) * 50_000 >>> boltons.statsutils.mean(data) 500000000.60091573 >>> sum(data) / len(data) 500000000.60091573
- MaxBarraclough 5y agoCan't we just use Kahan summation [0] then divide by n? edit Someone on the reddit thread [1] had the same thought. [0] https://en.wikipedia.org/wiki/Kahan_summation_algorithm https://en.wikipedia.org/wiki/Kahan_summation_algorithm [1] https://old.reddit.com/r/perl/comments/2zneat/how_you_average_numbers_matters/ https://old.reddit.com/r/perl/comments/2zneat/how_you_averag...
- m4r35n357 5y agoOr can we sum more accurately by taking the stable_mean and multiplying by n?
- SuchAnonMuchWow 5y agoKahan summation is better than naive summation, but not perfect. The stable_mean algorithm in the article is not perfect either. It just has more precision than the naive one, but it more or less has the same precision as twice-as-precise floats. >>> data = (1_000_000.1, 1.1, -1_000_000.1, -1.1) * 25_000 >>> sum(data) / len(data) == -2.3283153183228934e-16 >>> sum(sorted(data)) / len(data) == 2.7122441679239272e-11 >>> stable_mean(data) == -2.571225228716577e-13 >>> float(sum(numpy.array(data, dtype=numpy.float128)) / len(data)) == 2.2648549702353195e-19 Note that in python there is a math.fsum function which return a perfectly rounded summation of floats. There are multiple different techniques for that, but they are more complex than Kahan summation. >>> math.fsum(data) / len(data) == 0.0
- dataflow 5y agoSpecifically https://code.activestate.com/recipes/393090/ https://code.activestate.com/recipes/393090/
- Majromax 5y agoIn this case, it seems as if stable_mean is less accurate than the naive sum?
- SuchAnonMuchWow 5y ago
- amelius 5y agoDoes numpy.mean() do this?
- eesmith 5y agoSimple testing shows it does not use the naive sum. >>> import numpy as np >>> arr = np.array( (1_000_000_000.1, 1.1) * 50_000) >>> arr.mean() 500000000.6000003
- andi999 5y agoGive me an (real life) example where these errors matter. Any real data will be noisy anyway and relying on anything at the level of ppb (parts per billion) looks suspicious. Clarification: I do not say that treating floats careful is not important, and I know doing stuff like subtracting large numbers can lead to inaccurate leading digits. This is clear and part and parcel of any numerical work. But this is not what the article is about. I am claiming instead the article is looking for a problem where there is none (or let's say I fail to see) .
- SuchAnonMuchWow 5y agoRounding errors might be small, but catastrophic cancellation can make errors much larger than parts per billion. Maybe for your application it doesn't matter, but sometimes it does, and you need to be aware of these cases. You wouldn't want those kind of rounding errors in financial applications, in rocket science applications, etc ...
- andi999 5y agoAnd how does the mentioned averages procedure help?
- vardump 5y agoIf you're developing floating point code, you're in for a huge surprise, eventually. What happens for example when you subtract two close numbers with an error like this? There are a lot of mechanisms for small errors to become significant. Many of us have had their attitudes readjusted by a bug where these "insignificant" errors became anything but. We're often wrong in ways we can't imagine, so it is wise to remain humble and try to strive for correctness over our potentially misguided ideas how things work. For example, you'd better still have multiple layers of security even though we've already made some other layer "impenetrable". And still check for those "impossible" error cases before they create nearly undebuggable issues.
- andi999 5y agoYes, one has to careful about floating point, but honestly the article is not about that. Actually it is more about something where you do not have to be careful.
- jfrunyon 5y agoYAGNI. If you do, in fact, need to worry about floating point error, then there is not some magical "it will never happen if you have this much precision and use this neat algorithm". Heck, even fixed/decimal/rational numbers aren't perfect solutions, they just make it easier to tell if you're staying within your precision. Without infinite storage, you can't retain infinite precision, so the only real question is how much error you are willing to tolerate. And that depends very much on what you're doing - I can tolerate a lot more error in rolls on a MMO than in my bank account balance...
- johnbcoughlin 5y ago> the only real question is how much error you are willing to tolerate So, you are, in fact, gonna need an algorithm for summation that controls that error, if you're not willing to tolerate very much.
- jfrunyon 5y agoYes. Similarly, if you design something with total disregard to optimization, it will likely need optimization. YAGNI is a prediction, not a rule. Just like premature optimization, alternatives to floating point - or using longer floating point types - have real costs, which must be balanced against the needed precision.
- commandlinefan 5y agoWell, his final example suggests otherwise, but I kind of suspect that there must be something else going on there (that I'm too lazy to try to track down). Still, while fixing most floating point errors seem like rearranging the deck chairs on the titanic, if there's one algorithm that's 99.999% accurate and another that's only 99.9% accurate, might as well use and be aware of the more accurate one.
- jfrunyon 5y agoThe error in the 'Variance calculated using Statistics::Descriptive' is ~0.0002% of the mean (and of each value). As the OP says: "One saving grace of the real world is the fact that a given variable is unlikely to contain values with such an extreme range" My point is that 99.999% algorithm is going to have a cost over the 99.9% algorithm (more complex code, more resource usage), and depending on what you're doing, you may have to consider that cost.
- _0ffh 5y agoI suppose if you're not bothered with wasting lots of time, do the following for a precise sum: Select the two smallest numbers from the list, remove them and insert their sum. Repeat until only one number is left. [Ed. smaller meaning closer to zero]
- dolmen 5y agoThis is exactly what happens when you sort the number ascending and compute the sum starting from the smallest. Like in the article. Except this doesn't work if you must handle negative numbers. So what matters is not the smallest numbers but the number nearer from zero.
- MauranKilom 5y agoNo, that's not the same. Summing [4, 5, 6, 7] would yield the series (4, )9, 15, 22 with "sort then add", but what the previous poster described would first sum 4+5 (leaving the list [9, 6, 7]), then 6+7 (leaving the list [9, 13]) and finally arriving at 22.
- dolmen 5y agoOK. That algorithm could be improved (for accuracy, not speed) by replacing sums of the smallest number with multiplications as I suppose that multiplication of more than 3 identical numbers or more is loosing less precision than multiple additions. For example: [4, 3, 7, 7, 8] could be calculated as: (3+4) => [7, 7, 7, 8] (7*3) => [21, 8] (8+21) => 29
- metalliqaz 5y agoWow, Perl. Now, that's a name I've not heard in a long time.
- rjmunro 5y agoI still think the best way to deal with this stuff is to use interval arithmetic. Calculate the numbers with rounding set to round up, and separately to round down and compare the results. See PyInterval https://pyinterval.readthedocs.io/ https://pyinterval.readthedocs.io/ for an implementation. >>> from interval import interval >>> data = (interval(1_000_000_000.1), interval(1.1)) * 50_000 >>> sum(data) / len(data) interval([500000000.5971994, 500000000.601347]) The result is definitely between the 2 values. If you try another method to do the calculation, (e.g. sorting first) you will know which one was better without the answer having to be obvious as it is in this case. You might find one method is better at the lower bound, the other is better at the higher bound, which means you can narrow the value further than using either individually. It's a really powerful system that I think should be much better known.
- nestorD 5y agoYou should be okay with a sum but interval arithmetic can blow up in a way that is not correlated with the numerical accuracy of the algorithm. Newton's method, for example, is strongly convergent but produces diverging intervals.
- uyt 5y agoWhy does that happen? Is it because the errors are in opposite directions and will cancel each other out?
- nestorD 5y agoThat (and overall the unability to take correlation between error into account) and the fact that intervals propagated are upper bounds for some things like the division.
- rjmunro 5y agoYou can use other algorithms in these cases, and you can do things like converge your interval on a result by testing different values and seeing if the result is above the root or below the root, (or if it's within the error, so not possible to determine).
- jl2718 5y agotl;dr: sort first Better: recursive bucket accumulators using the exponent
- rhaps0dy 5y agoAccording to Radford Neal, it is possible to sum a sequence of n-bit floating point numbers exactly, i.e. rounding to the n-bit floating point value closest to the exact result. The algorithm uses 2n bits for the sum. https://radfordneal.wordpress.com/2015/05/21/exact-computation-of-sums-and-means/ https://radfordneal.wordpress.com/2015/05/21/exact-computati...
- asciimov 5y agoI have a special hatred of floating point numbers. I wish and advocate for a better standard. My issue is that I've come across errors due to floating point in weird places. One example, searching a list of strings for another string. If the sort function decides those strings look like numbers it may coerce them into doubles (ints being too small to hold them). Then sort wont match correctly due to the least significant bits being rounded off.
- Kranar 5y agoNot sure what this could possibly have to do with floating point numbers, nor am I aware of any function that does what you're describing. Is this some kind of JavaScript quirky behavior because nothing about floating point numbers has anything to do with coercing strings.
- TacticalCoder 5y agoI've definitely been bitten by similar (not identical) stuff. I've seen data serialized and sent over the wire by I-don't-know-what where several parsers can't correctly parse/serialize and end up with the same thing written back. Be it 32 bit float/64 bit float/string parsed as 32 bit float/string parsed as 64 bit float/string parsed as big decimals: I don't care, the fact is there are buggy serializers, buggy parsers, and even sometimes buggy parsers going into infinite loop (the JVM is famous for that: there was a bug where Double.parseDouble could be sent into an infinite loop by feeing it a carefully crafted string). The fact is: such kind of SNAFU does happen.
- spacechild1 5y agoI think the lesson is rather: don't use text based protocols for data serialization if you care about precision. Such things should never happen with a binary protocol.
- TacticalCoder 5y ago> If the sort function decides those strings look like numbers it may coerce them into doubles (ints being too small to hold them). You sometimes even see that with APIs inconsistently wrapping behind quotes or not (seemingly at random) some values. And then it gets trickier because you can have the parser correctly parse into, say, "big decimals" some of the values, while others may be automagically, imprecisely, parsed into floating-point numbers. The problem is of course made worse by people overall overusing decimal numbers where integers would be totally fine. One of the most beautiful floating-point bugs ever is this one (from about ten years ago as I recall it): language preference in browsers could be set as a decimal number between 0 and 1. Java (and other languages) had a broken floating-point implementation where parsing a certain floating-point number would lead to an infinite loop. So you could set one core of any webserver with that bug by simply requesting a page, with the language preference set to that floating-point number. Rinse and repeat a few times and all cores were spinning at 100%, effectively DDoS'ing the webserver.
- ball_of_lint 5y agoThat last example looks more like an outright bug than a floating point error.
- worewood 5y agoThis is useful for people working in critical applications. For the majority of people the 0.000001% difference doesn't matter.
- jart 5y agoThat's the wrong way to think about floating point. Double gives you ~15 digits of precision. If the answer is 500000000.600000 but your algorithm says 500000000.601239 then 25% of the bits you've got are kaput. It's also a question of scalability. Increase the author's example array size to 1 gigabyte and now 40% of your double is garbage data. Increase it to web scale, and you're toast. What's the solution? Divide and conquer. Not only will that let you farm out the summation across machines, but it also eliminates the linear growth in errors.
- worik 5y agoIf you are looking for such accuracy (where 50_000_000.00006 is different from 50_000_000.06) then you need a specialised library and do very careful coding. If you are averaging bazillions of numbers and they range from weeny to huge and you thing 50_000_000.000006 != 50_000_000.06 for your purposes then it is hard for me to see that you are doing anything sensible. This is, except in very specialised cases that I cannot imagine, a non issue. Just add them up and divide, (and check for overflow, that could ruin your day!)
- alex_smart 5y agoThe "Statistics::Descriptive" example (near the end of the article) has to be a bug in the relevant library, right? Even the naive floating implementation can't possibly get it this wrong for just 4 numbers.