[Mathe]: Wir konstruieren einen Raum auf Gauss-Kurven

[antisteo]

Hallo liebe ZFX+Developia-Community,

mich interessiert folgendes:

es gibt ja oft auch Werte (Double-Zahlen), dessen Wert wir nicht genau wissen. Ein Beispiel: Ich habe ein Foto und will daraus eine Punktwolke errechnen. Die Pixel geben zwar eine exakte Farbe an und anhand der Satelliten-Koordinaten bekomme ich auch X und Z heraus (und ein grobes Y), allerdings können immer wieder Abweichungen in den Werten auftreten, z.B.:

• Die Linse der Kamera weicht um einige tausendstel Grad ab
• Die Farbe ist durch Wolken verfälscht
• Ein zweites Bild, das darauf matchen soll, hat minimal andere Farbgebung
• und so weiter

Mir ist jetzt die Idee gekommen, anstatt mit "festen" Werten (double) mit Gauss-Kurven zu arbeiten, d.h. jede Werte-Angabe hat ein Zentrum (value) und eine Standardabweichung (deviation). Geschrieben werden diese z.B. als 4+/-0.2

Daraus bilden wir dann einen Raum, in dem man +,-,*,/ rechnen kann. Beispiel: 4+/-1 + 1 = 5+/-1. Durch das Addieren mit 1 erhöht sich die Abweichung nicht.

Eine Variable a=4 kann also nur den Wert 4 haben. Eine Variable a+/-0.1 kann hingegen JEDEN beliebigen Wert haben, allerdings sind diese Werte unterschiedlich wahrscheinlich. Der Vorteil: Man kann diese Werte in ein Gleichungssystem packen mit widersprüchlichen Aussagen in den Gleichungen und beliebig mehr Unbekannten als Werten (jede Standardabweichung ist ja quasi eine "unbenannte unbekannte") und das System lässt sich trotzdem lösen.

Jetzt die Preisfrage: Ich habe die 4 Grundrechenarten wie folgt definiert:

Code: Alles auswählen

// Gauss ist definiert als {value, deviation}
case "+":
   return Gauss{l.value + r.value, l.deviation + r.deviation}
case "-":
   return Gauss{l.value - r.value, l.deviation + r.deviation}
case "*":
   return Gauss{l.value * r.value, l.deviation * r.deviation + l.value * r.deviation + r.value * l.deviation}
case "/":
   return Gauss{l.value / r.value, l.deviation / r.value + l.value * r.deviation}

Folgende Aussagen funktionieren schon mit der Definition des Raums:

• 2 + 4 = 6
• 4+/-1 + 2 = 6+/-1
• 2 * 4+/-1 = 8+/-2
• 8+/-2 / 2 = 4+/-1

Ist das so plausibel? Fallen euch bessere Formeln ein?

[Krishty]

Ohne das im Detail gelesen zu haben: Mögliche Kollision mit Gaussian Splatting of Radiance Fields, das gerade sehr populär ist? https://repo-sam.inria.fr/fungraph/3d-g ... splatting/

Edit: Für Einsteiger

[Schrompf]

Die Formeln erinnern mich jetzt an die Fehlerrechnung, die bei den Praktika in meinem Studium immer gefordert waren. Und das isses ja auch, zumindest philosophisch: rechnen mit Unsicherheiten. Wenn Du das weiterbetreiben willst, kannst Du ja mal induktiv beweisen, dass Deine Version der Multiplikation aus der vielfach angewendeten Addition hervor geht.

Ansonsten seh ich da keine wirkliche Verwendung für.

[TomasRiker]

Das ist so ähnlich wie Intervallarithmetik, nur mit einer Gaußschen Verteilung statt begrenzter Intervalle. Ein Problem bei solchen Konstruktionen tritt auf, wenn derselbe Ausdruck mehrmals in einer Berechnung verwendet wird, denn dann kriegt man viel zu große Unsicherheiten. Z. B. x - x sollte eigentlich immer 0 ergeben, ohne jegliche Unsicherheit, aber du wirst stattdessen 0 mit einer größeren Unsicherheit rauskriegen.

[antisteo]

Das ist so ähnlich wie Intervallarithmetik, nur mit einer Gaußschen Verteilung statt begrenzter Intervalle. Ein Problem bei solchen Konstruktionen tritt auf, wenn derselbe Ausdruck mehrmals in einer Berechnung verwendet wird, denn dann kriegt man viel zu große Unsicherheiten. Z. B. x - x sollte eigentlich immer 0 ergeben, ohne jegliche Unsicherheit, aber du wirst stattdessen 0 mit einer größeren Unsicherheit rauskriegen.

Richtig. Alle in einer Formel vorkommenden Abweichungen müssen voneinander unabhängig sein, damit das Ergebnis aussagekräftig ist. An demselben Problem knabbere ich nämlich auch gerade, wo ich den Gleichungssolver bauen will. Da muss ich zuerst die Formel so umstellen, dass jede Variable nur einmal im Term vorkommt.

[antisteo]

Das ist so ähnlich wie Intervallarithmetik, nur mit einer Gaußschen Verteilung statt begrenzter Intervalle. Ein Problem bei solchen Konstruktionen tritt auf, wenn derselbe Ausdruck mehrmals in einer Berechnung verwendet wird, denn dann kriegt man viel zu große Unsicherheiten. Z. B. x - x sollte eigentlich immer 0 ergeben, ohne jegliche Unsicherheit, aber du wirst stattdessen 0 mit einer größeren Unsicherheit rauskriegen.

Danke für den heißen Tipp!

Die nutzen min() und max(), um die Abweichung einzugrenzen, ich hatte mir das nur nicht getraut

[Jonathan]

Man muss auch bedenken, dass Abweichungen als Gausskurve zu beschreiben, in gewisser Weise ein naiver (im Sinne von "erst-bester") Ansatz ist. Viele Fehler sind zwar Normalverteilt (deswegen "erst-bester" und nicht "erst-beliebiger"), aber das muss nicht so bleiben wenn man damit weiter rechnet.

Beispiel: Du dividierst 20 durch den unsicheren Wert 10 mit einer Standardabweichung von 5. Der Erwartungswert ist, dass 20/10=2 raus kommt, es kann aber gleichwahrscheinlich auch 20/5=4 oder 20/15=1.33 herauskommen. Die beiden Zahlen sind nun aber unterschiedlich weit vom Erwartungswert entfernt, die Verteilung aller möglichen Werte ist also nicht mehr symmetrisch, auch wenn der Fehler des Messwertes 10 ursprünglich perfekt symmetrisch und Normalverteilt war.
Die wahre neue Verteilung kannst du jetzt unmöglich mit 2 Werten (Mittelwert und Abweichung) beschreiben. Du kannst die Gesamtverteilung zwar erneut mit einem Gauss approximieren, aber kriegst dadurch zusätzliche Unsicherheiten rein die du nicht vermeiden kannst. Zwangsläufig wird der Fehler überschätzt, unterschätzt oder beides gleichzeitig, aber nie korrekt geschätzt.

Das ganze wird entsprechend schwieriger, wenn du auch Wurzeln und Potenzen und Dinge wie sinus-Funktionen zulassen willst. Stell dir einfach mal vor du berechnest die Wurzel einer Zahl die du nur sehr ungenau kennst. Sagen wir, das Ergebnis ist "2, mit Standardabweichung 4". Aber das Ergebnis einer Wurzel ist strikt positiv, -2 wäre also immer unmöglich, und kann in einigen Formeln Probleme bereiten. Du kannst das schlicht nicht mit einem Gauss abbilden.

Im Allgemeinen ist es denke ich so, dass deine Beschreibung des Fehler so kompliziert ist, wie die ursprünglichen Verteilungen und die Berechnungen die du damit anstellst. Sprich, was du versuchst mit 2 Zahlen auszudrücken, könnte in Wirklichkeit hunderte Zahlen erfordern. Deine Modellierung von Ungenauigkeiten führt also selber zusätzliche Ungenauigkeiten ein.

[antisteo]

Moin,

erst mal Danke für die Eingaben. Ich weiß genau, dass die Normalverteilung nur eine Näherung ist. Allerdings interessieren mich hier nur lineare Gleichungssysteme und bei denen ist ja alles so schön linear. Im Prinzip will ich Kamera-Positionen von Video-Frames erraten. Da sind ein paar Matrix-Multiplikationen mit einer unbekannten Matrix und ein paar Millionen lineare Lichtstrahlen, die aus der Kamera kommen und in einer Pixelfarbe enden.

Der Trick mit dem max() aus deviation hat Wunder bewirkt! Die Ergebnisse meiner Annäherungen werden auf einmal viel genauer, die Erwartungswerte treffen sogar genau ins Schwarze. Wenn man sich das grafisch vorstellt, ergibt das auch Sinn. Legt man zwei Gauss-Kurven mit unterschiedlicher Standardabweichung übereinander, wird die kleinere von der größeren Kurve "verschluckt".

Was auch noch hinzukommt: Weiß man am Ende jeden Pixel inkl. seiner Position und dessen Standardabweichung, hat man auch noch die Punktwolken-Größe (Punktgröße = Standardabweichung)

[antisteo]

Habemus fertig

Folgende Sammlung an Gleichungssystemen:

Code: Alles auswählen

a = 4+/-0.5
b = a + 2
c = a + b
#c = 2 * a + 1

f = 2*x+1+/-0.5 # first equation
f = -4*(x-1) - 1 # second equation
# solution is: x = 1/3, f = 1.666
#x = 0.5+/-100

w / 5 - 1

ergibt folgenden Output des Solvers:

Code: Alles auswählen

starting expansion
2 * x + 1+/-0.5 = -4 * (x - 1) - 1
  x = 0.3333333333333333+/-0.08333333333333333
  0 = 6 * (2 * x + 1+/-0.5 - (-4 * (x - 1) - 1)) + (2 * x + 1+/-0.5 - (-4 * (x - 1) - 1)) * 6
-> x = 0.3333333333333333+/-0.041666666666666664

2 / (f - 1+/-0.5) = 1 + -4 / (1 + f)
  0 = (-2 / (f - 1+/-0.5 * (f - 1+/-0.5)) - 4 / (1 + f * (1 + f))) * (2 / (f - 1+/-0.5) - (1 + -4 / (1 + f))) + (2 / (f - 1+/-0.5) - (1 + -4 / (1 + f))) * (-2 / (f - 1+/-0.5 * (f - 1+/-0.5)) - 4 / (1 + f * (1 + f)))
4+/-0.5 = b - 2
  b = 6+/-0.5
  0 = -1 * (4+/-0.5 - (b - 2)) + (4+/-0.5 - (b - 2)) * -1
-> b = 6+/-0.25

4+/-0.5 = c - b
  c = b + 4+/-0.5
  0 = -1 * (4+/-0.5 - (c - b)) + (4+/-0.5 - (c - b)) * -1
-> c = -(-1 * (4+/-0.5 - -b) + (4+/-0.5 - -b) * -1) / 2

  b = c - 4+/-0.5
  0 = 4+/-0.5 - (c - b) + (4+/-0.5 - (c - b))
-> b = -(4+/-0.5 - c + (4+/-0.5 - c)) / 2

b - 2 = c - b
  c = b + (b - 2)
  0 = -1 * (b - 2 - (c - b)) + (b - 2 - (c - b)) * -1
-> c = -(-1 * (b - 2 - -b) + (b - 2 - -b) * -1) / 2

  b = -(-2 - c) / 2
  0 = 2 * (b - 2 - (c - b)) + (b - 2 - (c - b)) * 2
-> b = -(2 * (-2 - c) + (-2 - c) * 2) / 8

a + 2 = c - a
  a = -(2 - c) / 2
  0 = 2 * (a + 2 - (c - a)) + (a + 2 - (c - a)) * 2
-> a = -(2 * (2 - c) + (2 - c) * 2) / 8

  c = a + (a + 2)
  0 = -1 * (a + 2 - (c - a)) + (a + 2 - (c - a)) * -1
-> c = -(-1 * (a + 2 - -a) + (a + 2 - -a) * -1) / 2

a - 4+/-0.5
b - (a + 2)
c - (a + b)
f - (2 * x + 1+/-0.5)
f - (-4 * (x - 1) - 1)
w / 5 - 1

b=a + 2
b=c - a
b=6+/-0.5
b=6+/-0.25
b=c - 4+/-0.5
b=-(4+/-0.5 - c + (4+/-0.5 - c)) / 2
b=-(-2 - c) / 2
b=-(2 * (-2 - c) + (-2 - c) * 2) / 8
---
b=5.999999996812187+/-4.440892098595651e-16

c=a + b
c=b + 4+/-0.5
c=-(-1 * (4+/-0.5 - -b) + (4+/-0.5 - -b) * -1) / 2
c=b + (b - 2)
c=-(-1 * (b - 2 - -b) + (b - 2 - -b) * -1) / 2
c=a + (a + 2)
c=-(-1 * (a + 2 - -a) + (a + 2 - -a) * -1) / 2
---
c=9.999999993624373+/-3.8010272800181775e-26

f=2 * x + 1+/-0.5
f=-4 * (x - 1) - 1
---
f=1.666666666666667+/-0.16666666666666666

x=2 / (f - 1+/-0.5)
x=1 + -4 / (1 + f)
x=0.3333333333333333+/-0.08333333333333333
x=0.3333333333333333+/-0.041666666666666664
---
x=0.33333333333333326+/-0.04166666666666672

w=5
---
w=5

a=4+/-0.5
a=b - 2
a=c - b
a=-(2 - c) / 2
a=-(2 * (2 - c) + (2 - c) * 2) / 8
---
a=3.999999996812186+/-4.440892098595651e-16

in anderen Worten: Ich kann jetzt beliebige lineare Gleichungssysteme in den Solver packen und er berechnet mir trotz bekannter Ungenauigkeiten in den Messdaten die Parameterbelegung mit der höchsten Wahrscheinlichkeit.

Als nächstes muss ich anfangen, die umgestellten Formeln in ein Pytorch-Skript zu generieren und massenhaft Bild- und Videodaten da reinzuschicken und erhalte hoffentlich hübsche Punktwolken von der Welt da draußen.

Der Solver kann unter anderem:

• Terme im Gauss-Körper ausrechnen
• Terme nach Variablen umstellen
• Ableitungen nach einer Variablen bilden
• Polynome lösen
• Terme mit der Nullpunkt-Methode lösen (Term quadrieren, Ableiten und die Nullstelle in der Ableitung ist die Lösung)

[TomasRiker]

Wo du "quadrieren" erwähnst: Wie bildest du das ab, dass ein Ausdruck nach dem Quadrieren nicht negativ sein kann? Das ist ja dann keine Normalverteilung mehr.
Nur das Produkt zweier unabhängiger Normalverteilungen ist wieder eine Normalverteilung. Quadrieren bedeutet ja mit sich selbst multiplizieren, also nicht unabhängig.

[Jonathan]

in anderen Worten: Ich kann jetzt beliebige lineare Gleichungssysteme in den Solver packen und er berechnet mir trotz bekannter Ungenauigkeiten in den Messdaten die Parameterbelegung mit der höchsten Wahrscheinlichkeit.

Warte, es ging die ganze Zeit nur darum, die Lösung mit der höchsten Wahrscheinlichkeit zu finden? Unter der Annahme, dass alle Fehler immer Normalverteilt sind? Dann ist der Erwartungswert ja immer die Mitte, d.h. du kannst einfach alle Fehlerabschätzungen weglassen und mit den Werten rechnen, als seien sie exakt, und es wird immer das selbe rauskommen. So scheint mir jetzt aber, dass du viel Aufwand treibst, den du dann am Ende wegschmeißt. Sinnvoll wird die ganze Mühe erst, wenn du irgendwann mal kompliziertere Verteilungen benutzt...

[antisteo]

in anderen Worten: Ich kann jetzt beliebige lineare Gleichungssysteme in den Solver packen und er berechnet mir trotz bekannter Ungenauigkeiten in den Messdaten die Parameterbelegung mit der höchsten Wahrscheinlichkeit.

Warte, es ging die ganze Zeit nur darum, die Lösung mit der höchsten Wahrscheinlichkeit zu finden? Unter der Annahme, dass alle Fehler immer Normalverteilt sind? Dann ist der Erwartungswert ja immer die Mitte

nein, ist er nicht. Die vielen Werte-Vorschläge werden nach der Wahrscheinlichkeit noch mal gewichtet.

Für Werte \(x_i\) als Wert und \(\sigma{}_i\) als Standardabweichung:

\(x = \Sigma{}_i (x_i \Pi{}_{j, i \neq j} {G(x_i, x_j, \sigma{}_j)}) / \Sigma{}_i \Pi{}_{j, i \neq j} {G(x_i, x_j, \sigma{}_j)}\)
\(\sigma = \min{}_i (\sigma{}_i + |x - x_i|)\)

, d.h. du kannst einfach alle Fehlerabschätzungen weglassen und mit den Werten rechnen, als seien sie exakt, und es wird immer das selbe rauskommen. So scheint mir jetzt aber, dass du viel Aufwand treibst, den du dann am Ende wegschmeißt. Sinnvoll wird die ganze Mühe erst, wenn du irgendwann mal kompliziertere Verteilungen benutzt...

nope, da ich brauche die Standardabweichung tatsächlich. Zum einen weiß ich mit der Standardabweichung auch gleichzeitig die Größe des Farbkleckses auf dem Bildschirm (wenn man jetzt Gauss Splats benutzt). Zum anderen gibt es immer mehrere Lösungsvorschläge je Wert (zwei Pixel aus zwei verschiedenen Bildern zeigen auf dieselbe Weltkoordinate) - mit der Standardabweichung weiß ich, welche der Lösungen am richtigsten ist. Im Prinzip sind bei einem Gauss mit \(\sigma > 0\) alle Lösungsvorschläge richtig - aber eben entsprechend unwahrscheinlich je weiter sie von der Norm entfernt sind.