Interview-Frage [asking for a friend]: Mittelwert von Array von Double

[Alexander Kornrumpf]

Moin,

ich habe ungefähr soviel Numerik gemacht, dass ich weiß, dass es nicht ganz so einfach ist, ohne (immer) zu wissen, wie es richtig geht.

Heute habe ich gelernt, dass ich die Antwort auf eine sehr einfache Frage nicht weiß: Wie berechnet man korrekt den Mittelwert ("Durchschnitt") eines Array von Double?

Es ist eine Interview-Frage, aber keine die mir gestellt wurde, es macht mich einfach nur neugierig und das hier ist ja die passende Zielgruppe dafür.

Naiv summieren und dann teilen funktioniert nicht, weil die laufende Summe irgendwann groß gegen die Summanden wird und dann bekommst du Genauigkeitsprobleme. Eine mögliche Lösung dafür ist https://en.wikipedia.org/wiki/Kahan_summation_algorithm, was in meiner Masterarbeit vorkam, deswegen werde ich das nie vergessen, ich hatte es hier im Forum auch schonmal erwähnt.

Das hilft aber natürlich auch nur begrenzt weiter. Der Interviewer, ich war ja nicht dabei aber so wurde es mir berichtet, hat speziell danach gefragt wie man verhindert, dass der Wertebereich von Double überschritten wird. Ich kann jetzt nicht beschwören, dass er mit Wertebereich nicht doch vielleicht Genauigkeit meinte, als Interviewee hätte ich das nachgefragt. Lass uns aber spaßeshalber annehmen wir haben ein sehr langes Array und die Summe läge tatsächlich über max_double. [Ich persönlich glaube immer noch dass Genauigkeit in der Praxis lange vorher ein Problem werden wird, vielleicht war es auch eine Fangfrage, oder der Interviewer hat es selbst nicht so ganz kapiert, keine Ahnung]

Jetzt sind wir also "schlau" und teilen erst und summieren dann, dann ist die Summe kleiner, sollte man denken. Aber um das zu machen muss man eine Alltagszahl durch eine sehr große Zahl (die Länge des Array) teilen und das ist auch eine schlechte Idee.

Also machen wir was?

Ich habe den Fehler gemacht das zu googlen und die Stackoverflow-Antworten, die ich gesampled habe waren alle von der Art dass sie das Problem gar nicht kapiert haben (oder kapiere ich es nicht?) aber std::accumulate kannten und stolz darauf waren. Ist natürlich auch schwer das im Public Internet zu fragen weil es ja sofort als Interview / Homework Toy-Problem zu erkennen ist, aber ich glaube irgendwie nicht, dass die Leute, die es als solches abtun, wirklich echt die Lösung kennen.

Oder stelle ich mich gerade dumm an und merke es nicht?

[Alexander Kornrumpf]

Ich habe es jetzt nicht bewiesen, aber Intuition, vage Erinnerung und ein bisschen Rumspielen mit Excel sagen beim zweiten Nachdenken, dass der Mittelwert linear ist (Summe der Mittelwerte ist der Mittelwert der Summe), man kann also das sehr große Array in kleinere Teilpakete aufteilen, die klein genug sind, dass die Summe noch in double passt, dann deren Mittelwerte ausrechnen, und dann den Mittelwert der Mittelwerte. Die Pakete müssen halt gleich groß sein. Zur Not macht man das rekursiv. Ich bin ehrlich gesagt gerade zu faul das mal auf einer Serviette durchzu-x-en, riskiere lieber mich zu blamieren :)

[starcow]

Vielleicht verstehe ich deine Frage jetzt nicht ganz richtig...
Doch du könntest ja einfach die Differenz zwischen dem nächsten Element und deinem aktuellen Mittelwert berechnen, diese entsprechend gewichten (wenns das dritte Element wäre, würde die Differenz durch 3 geteilt) und dann addierst du diesen Wert zu deinem ursprünglichen Mittelwert.
Bei einem 64Bit double hast du ja 11 Bit für den Exponent, was bedeuten würde, du hättest konkret für den Exponent -1022 (also 2^-1022) nach IEEE 754 (denormalisierte Zahlen mal weggelassen). Selbst wenn dein Array mehrere Milliarden Einträge hätte, dürftest du deswegen eigentlich noch lange nicht in Bedrängnis geraten (punkto Genauigkeit).

Edit:
Oder wie du ja selbst vorschlägst, machst du zuerst Teilpakete und berechnest dann jeweils deren Mittelwerte.

Edit 2:
Wenn ich etwas länger darüber nachdenke:
Die Frage nach dem Overflow ist relativ trivial (in den Interviews wirds wohl darum gehen) - "deine" Frage, nach der Genauigkeit ist wohl etwas anspruchsvoller!

[Schrompf]

Gibt doch so ein Rolling Average, wo man jeden Wert mit sinkendem Gewicht einberechnet und den bisherigen Durchschnitt anteilig. Tödlich für die Genauigkeit, aber dafür überlaufsicher.

Code: Alles auswählen

double average = 0;
double count = 0; // zu faul, im Beispiel nen Integer zu nehmen.
for( double d : values ) {
  ++count;
  average = (average * (1.0 - 1.0 / count)) + (d * 1.0 / count);
}

[starcow]

Müsste man dann eigentlich, bei so vielen Divisionen (mehrere Milliarden Elemente), konsequenter nicht auch noch die Frage nach der Performance aufwerfen? Gegenüber simplen Additionen müsste dies ja irgendwann mal mal ins Gewicht fallen, oder?

[NytroX]

Ich musste grad auf Wikipedia schauen, was du mit "Mittelwert" meinst.
Daran erkennt man den verkorksten IT-ler...

Du meinst offenbar den Durchschnitt, also das arithmetische Mittel, und nicht den Median oder das gewichtete oder geometrische oder quadratische Mittel :-)

Kann man nicht einfach ein binary128 nehmen anstatt ein binary64 (a.k.a. double)? Je nach Sprache und Compiler halt...
(Was ist eigentlich der name dafür? "quadle" ?)

Oder man könnte 2 doubles nehmen um die doppelte Präzision zu erhalten, indem man den Präzisionsfehler fortlaufend korrigiert (double-double Arithmetik):
https://www.youtube.com/watch?v=6OuqnaHHUG8 bzw.
https://www.youtube.com/watch?v=5IL1LJ5noww

Aber vielleicht reicht auch einfach double (oder long double, also 80bit), ich meine die Werte sind eh nicht präzise, sonst hätte man ja double nicht als Typ genommen...

[Krishty]

Ich habe es jetzt nicht bewiesen, aber Intuition, vage Erinnerung und ein bisschen Rumspielen mit Excel sagen beim zweiten Nachdenken, dass der Mittelwert linear ist (Summe der Mittelwerte ist der Mittelwert der Summe), man kann also das sehr große Array in kleinere Teilpakete aufteilen, die klein genug sind, dass die Summe noch in double passt, dann deren Mittelwerte ausrechnen, und dann den Mittelwert der Mittelwerte. Die Pakete müssen halt gleich groß sein. Zur Not macht man das rekursiv. Ich bin ehrlich gesagt gerade zu faul das mal auf einer Serviette durchzu-x-en, riskiere lieber mich zu blamieren :)

Noch zwei Gedanken:

1. Mit Gleitkommazahlen ist Division durch / Multiplikation mit Zweiterpotenz verlustfrei (von Über- und Unterlauf abgesehen). Die zur Not rekursive Methode könnte, sofern sie auf Paketen der Größe von Zweiterpotenzen arbeitet, also wirklich ziemlich genau sein. Bei einem Binärbaum für das Mittel von 1024 Zahlen würdest du die Wurzel mit maximal 10 ULP Abweichung erreichen(?)

2. Wenn du doch noch Über- und Unterlauf sinnvoll ausschließen möchtest, ist das Addieren von zwei Zahlen auf einmal wohl wirklich das Höchste Praktikable. Schließlich versagt (a + b) / 2 bei 2× DBL_MAX und a/2 + b/2 bei 2× DBL_MIN. Ich tendiere also wirklich zu was Binärbaum-Ähnlichem.

ich meine die Werte sind eh nicht präzise, sonst hätte man ja double nicht als Typ genommen...

Finde ich extrem verkehrt ;)

[Alexander Kornrumpf]

Danke für die vielen Antworten. Ich versuche mal ein kohärentes Ganzes zu antworten, statt jetzt alle einzeln zu zitieren.

Ich glaube Antworten, die darauf abzielen, dass man in der Praxis nicht in diese Situation käme oder einen größeren Datentyp nehmen könnte, sind vielleicht richtig, aber wir erinnern uns, dass das mit Absicht ein Toy-Problem ist das irgendwas testen will, vermutlich wie gut sich der Proband mit Fließkommazahlen auskennt. Es ist aber auch nicht völlig praxisfern, an die Genauigkeitsgrenzen kommt man nach meiner Erfahrung sehr schnell. Ich fage mich halt ob die Frage nach dem Wertebereich eine Fangfrage ist um davon abzulenken, dass Genauigkeit das eigentliche Problem ist, oder ob der Interviewer so weit gar nicht gedacht hat. Wir werden es nie erfahren.

Ich denke auch, dass in dem Fantasieuniversum in dem diese Aufgabe spielt, die einzige Hardwarelimitation ist dass wir mit 64bit Fließkomma auskommen müssen, aber von denen können wir beliebig viele haben und sie auch beliebig schnell verrechnen. So würde ich das zumindest als Interviewer framen (da ich jetzt mit meinem Klarnamen damit im Internet stehe, kann ich die Frage sowieso nie im Interview stellen, es gibt Spaßvögel die googlen den Interviewer vorher).

Der Hinweis von Krishty, dass man den Fall mit dem großen Array mit meinem eigenen "Trick" auf den Fall mit zwei Zahlen reduzieren kann macht einiges klarer.

[Lord Delvin]

Naja als Interviewfrage finde ich das gar nicht so schlecht, weil du schon sehen kannst, wie Leute über so Probleme denken oder was sie so erlebt haben. Meine Reaktion wäre btw. völlig anders gewesen, deswegen noch ein paar kurze Kommentare:

Vollkommen Praxisfern ist das nicht; kann z.B. transformation einer Metrik in einen Mittelwert für ein Dashboard sein. Die können mit 100MBit/s reinkommen.

Meine Antwort wäre "xs.fold(0)(_+_)" gewesen und dann hätte ich diskutiert, was genau das Problem damit ist; einfach weil meine Erfahrung ist, dass man aufhören kann, wenn in xs die *richtigen* doubles drinstehen; meine eigene Kritik an meiner Lösung wären Inf und NaN. Das würden eure Vorschläge vermutlich auch allesamt erstmal nicht richtig machen. *Richtig* ist da ohne Gegenfrage auch nicht drin.

Mir ist unklar, wie sich der numerische Fehler verhält wenn man die Foldlösung auf einer 80bit FPU hat; wenn das in der Praxis ein Problem ist würde ich es an einen Numeriker delegieren; im Internet zu suchen kann bei so speziellen Problem zu extrem gruseligen Ergebnissen suchen (habe in den letzten Wochen was im Detail über base64 wissen müssen; ist komplett zum Kotzen was man da im Internet zu findet)

Binärbaum fände ich seltsam, würde eher Richtung 2^{14-16} schauen; erfordert Gegenfragen und ggf. in der Realität Experimente; jenachdem was da das Problem ist würde ich mit Threads anfangen, btw. bevor ich Assembler verwenden würde; das liegt an meinen Fähigkeiten und meinem derzeitigen Arbeitsumfeld; hier muss einem natürlich klar sein, dass das Array ggf. nicht PowerOfTwo Länge hat.

Hätte auf die Gegenfrage gestellt, warum man nicht fixpunkt verwendet, wenn das ein Problem ist. Da muss man sich aber vorher überlegen auf was man sich genau beworben hat; Fin bestimmt ein guter Punkt, Games würde ich das nicht machen :)

[Krishty]

Mir ist unklar, wie sich der numerische Fehler verhält wenn man die Foldlösung auf einer 80bit FPU hat

Welche Sprache verwendet denn heute noch die 80-Bit-FPU?

[udok]

Wenn da nur ein einziges INV oder NaN drinnen ist, hat sich die Frage erübrigt... man muss also die Statistik ungefähr kennen, und den Anwendungszweck.
Wenn Geschwindigkeit wichtig ist, dann braucht man einen größeren Hardware Akkumulator, man könnte auch alle 1000 Elemente die Summe runterskalieren.
Wenn die Daten im 32 bit Float Format vorliegen, nimmt man double für die Summe, früher wäre Long Double für Daten im Double Format eine Option gewesen, wobei es den x87 ja noch immer gibt.
Wenn die Zahlen stark unterschiedliche sind, und Geschwindigkeit nicht wichtig ist, könnte man das Array vorher sortieren und mit den kleinen Werten anfangen.

[Lord Delvin]

Mir ist unklar, wie sich der numerische Fehler verhält wenn man die Foldlösung auf einer 80bit FPU hat

Welche Sprache verwendet denn heute noch die 80-Bit-FPU?

Ok, bin vielleicht jetzt zu weit weg davon, aber ich hatte irgendwie im Kopf, dass normale float-Arithmetik auf Intel CPUs in 80-Bit-FPUs verschoben wird, wenn sie nicht Schleifenoptimiert wurde. Bei der Gegenfrage würde ich jetzt erwarten, dass sich da was geändert hat :)

[udok]

Heute wird die x87 mit den 80 Bit Floats nur mehr in Legacy Apps verwendet. Jeder aktuelle Compiler erzeugt SSE2 Floating Point Code , und da gibt es nur 32 und 64 Bit Floating Point Formate.

[Alexander Kornrumpf]

Naja als Interviewfrage finde ich das gar nicht so schlecht, weil du schon sehen kannst, wie Leute über so Probleme denken oder was sie so erlebt haben.

Sehe ich auch so.

Meine Antwort wäre "xs.fold(0)(_+_)" gewesen und dann hätte ich diskutiert, was genau das Problem damit ist; einfach weil meine Erfahrung ist, dass man aufhören kann, wenn in xs die *richtigen* doubles drinstehen;

OK, was ist xs?

meine eigene Kritik an meiner Lösung wären Inf und NaN. Das würden eure Vorschläge vermutlich auch allesamt erstmal nicht richtig machen.

Ich glaube das ist wirklich out of scope für die Frage, würde aber in einem Interview vermutlich gerne gesehen.

im Internet zu suchen kann bei so speziellen Problem zu extrem gruseligen Ergebnissen suchen (habe in den letzten Wochen was im Detail über base64 wissen müssen; ist komplett zum Kotzen was man da im Internet zu findet)

Wie oben geschrieben kann man selbst zu der allereinfachsten Version des Problems nur std::accumulate "Experten" auf Stack finden, davon aber reichlich. Wenn man nicht schon weiß was das Problem ist copy und pasted man das halt in die Codebasis und gibt es dann an den nächsten Frager weiter.

Binärbaum fände ich seltsam, würde eher Richtung 2^{14-16} schauen; erfordert Gegenfragen und ggf. in der Realität Experimente; jenachdem was da das Problem ist würde ich mit Threads anfangen, btw. bevor ich Assembler verwenden würde; das liegt an meinen Fähigkeiten und meinem derzeitigen Arbeitsumfeld; hier muss einem natürlich klar sein, dass das Array ggf. nicht PowerOfTwo Länge hat.

Schlimmstenfalls, wenn die Länge des Arrays eine Primzahl ist, muss man bei jedem Teiler den Rand irgendwie behandeln. Ich führe selbst viele Interviews, die guten Kandidaten sind die, die solche Probleme erkennen und benennen können, aber ich würde raten, immer mit der einfachsten Version des Problems anzufangen und zu schauen, wohin der Kandidat rennt. Das relevante an der Krishty-Antwort für die fiktive Interview Situation ist mMn dass sie den komplizierten Fall auf den einfachst-möglichen zurückführt, nicht dass man sie ohne Stolpersteine implementieren kann. [Sorry für den "I'm a hiring manager AMA"-Modus, vielleicht hilft es irgendeinem Leser]

Hätte auf die Gegenfrage gestellt, warum man nicht fixpunkt verwendet, wenn das ein Problem ist. Da muss man sich aber vorher überlegen auf was man sich genau beworben hat; Fin bestimmt ein guter Punkt, Games würde ich das nicht machen :)

Du bist vielleicht auch nicht der ganz typische Bewerber.