Mathematische Rätsel mit bedingten Wahrscheinlichkeiten knackt der Fachmann mit der Bayes-Formel oder auch mit diskreten Verteilungen, erzeugt von kurzen Perl-Skripts.
Mit dem in Statistikerkreisen bekannten Ziegenproblem [6] und seiner verblüffenden Lösung kann der Rätselfreund reihenweise mathematische Laien hereinlegen, aber auch so mancher Fachmann hat sich daran schon damit blamiert, weil er leidenschaftlich gegen die richtige Lösung gewettert hat. Wer hätte gedacht, dass sich Wahrscheinlichkeiten in offensichtlich festgezimmerten Fernsehstudios dramatisch ändern, nur weil der Moderator eine weitere Tür ohne Preis öffnet?
Das menschliche Gehirn scheint sich allgemein schwer zu tun mit sogenannten bedingten Wahrscheinlichkeiten, die zufällige Ereignisse mit Vorbedingungen beschreiben. So auch bei der vor 15 Jahren im Perl-Snapshot ([4]) schon einmal vorgestellten Aufgabe: "Gegeben ist ein Zylinder, in dem drei Karten liegen: Die erste ist vorne und hinten schwarz, die zweite vorne und hinten rot, und die dritte auf der einen Seite schwarz und auf der anderen rot. Eine Karte wird gezogen. Man sieht die Vorderseite, und die ist rot. Wie hoch ist die Wahrscheinlichkeit daß auch die Rückseite der gezogenen Karte rot ist?"
Abbildung 1: Der Satz von Bayes (von der Wikipedia) |
Die meisten Leute antworten darauf mit "50%", denn scheinbar gibt es zwei gleichwahrscheinliche Fälle, die rot-schwarze und die rot-rote Karte, und einmal ist die Kehrseite schwarz und einmal rot. Steckt der Rätselfreund dem erstaunten Zuhörer aber dann, dass die korrekte Lösung "66%" lautet, macht sich oft Entrüstung breit. Der Kniff liegt an den Vorbedingungen des Experiments: Die rot-rote Karte wird einfach doppelt so häufig gezogen wie die rot-schwarze, weil das Experiment als Vorbedingung "eine Seite ist rot" stellt.
Abbildung 2: Aus einem Zylinder mit einer schwarz/roten, einer schwarz/schwarzen und einer rot/roten Karte zieht ein Proband Karten. |
Schon vor 250 Jahren hat sich der Mathematiker Thomas Bayes daran gemacht, solche gefühlsmäßig schwer zu verstehenden Probleme in eine Formel zu fassen, an der sich nicht rütteln lässt. Der "Satz von Bayes" ([5]) beschreibt in seiner "diachronischen Interpretierung", wie sich die Wahrscheinlichkeiten von Hypothesen im Lauf der Zeit verändern, falls in einem Experiment neue Daten auftauchen. Die Formel
P(H|D) = P(H) * P(D|H) / P(D)
definiert die Wahrscheinlichkeit einer Hypothese H (zum Beispiel "ich habe die rot-rote Karte gezogen") in Abhängigkeit von neu auftauchenden Daten D ("die Vorderseite der gezogenen Karte ist rot"). Auf der rechten Seite ist P(H) dabei die Ausgangswahrscheinlichkeit der Hypothese, noch bevor neue Daten vorliegen. Multipliziert wird sie mit P(D|H), also der Wahrscheinlichkeit, dass die Daten unter der Hypothese vorliegen. Wie wahrscheinlich ist es, dass der Kartenfreund auf eine rote Vorderseite starrt, falls er die rot-rote Karte gezogen hat? Und im Nenner der rechten Seite steht schließlich mit P(D) die Wahrscheinlichkeit der vorliegenden Daten, unabhängig von irgendeiner Hypothese. Wie wahrscheinlich ist es, dass der Proband nach dem Ziehen irgendeiner Karte auf eine rote Vorderseite starrt?
Alle drei möglichen Hypothesen (zufällig gezogene Karten) sind in diesem Experiment offensichtlich gleich wahrscheinlich. Also ist P(H), also die Wahrscheinlichkeit, dass die rot-rote Karte gezogen wird, gleich 1/3, denn sie wird in einem Drittel aller Fälle gezogen. Und bei Annahme dieser Hypothese ist die Kehrseite der gezogenen rot-roten Karte trivialerweise in exakt 100% aller Fälle ebenfalls rot, also ist P(D|H) = 1. Unabhängig von der Hypothese einer gezogenen Karte ist die Wahrscheinlichkeit, nach dem Ziehen auf eine rote Kartenoberseite zu blicken 50%, also ist P(D) = 1/2. Schließlich liegen im Zylinder 6 Kartenhälften, von denen drei rot und drei schwarz sind. In die Bayes-Formel eingesetzt ergibt sich für P(H|D) entsprechend 1/3 * 1 / (1/2) = 2/3. Die Bayes-Formel sagt also mit 66% das im Experiment nachweisbar korrekte Ergebnis voraus, und schlägt damit der falsch liegenden Intuition ein Schnippchen.
Das etwa ein halbes Jahr alte Buch "Think Bayes" [2] beschreibt nun zusätzlich zur Bayes-Formel ein numerisches Verfahren, das sich zum Experimentieren leicht in handliche Skripts fassen lässt. Dazu legt man die Hypothesen in einer statistischen Verteilung ab, speichert sie also in einem Python-Dictionary oder Perl-Hash mit ihren Ausgangswahrscheinlichkeiten ab, multipliziert dann die Einzelwerte mit den bedingten Wahrscheinlichkeiten eintreffender Daten unter Annahme der jeweiligen Hypothese, normalisiert dann die Werte für alle Hypothesen, und erhält die Wahrscheinlichkeit der gesuchten Hypothese im Lichte neu eingetroffener Daten. Dieses Verfahren, und speziell die abschließende Normalisierung funktioniert allerdings nur, falls:
1) höchstens eine der Hypothesen wahr ist und
2) es keine weiteren Möglichkeiten gibt, also mindestens eine der Hypothesen wahr ist
Listing 1 zeigt die Implementierung des Skripts cards
, das am Ende
die gesuchte Lösung "0.666.." ausgibt. Es nutzt das in Listing 2
vorgestellte Modul Distrib.pm
, das die Hilfsfunktionen zur Berechnung
der statistischen Verteilung bereitstellt. Die Wahrscheinlichkeiten
für die Hypothesen der Ziehung verschiedener Karten
("BB", "BR", "RR" für black-black, black-red, red-red)
füttern die Zeilen 7-9 mit jeweils dem Wert 1/3 in das Modul.
Die hinzugekommenen Daten, also die Oberseite der tatsächlich
gezogenen Karte, geben die Zeilen 11-13 für jede Hypothese vor und
multiplizieren die eingespeisten Hypothesenwerte damit. So beträgt die
Wahrscheinlichkeit, dass sich dem Proband eine rote Oberseite präsentiert
bei der schwarz-schwarzen Karte gleich Null. Bei der schwarz-roten Karte
beträgt sie 50%, je nachdem wie herum die Karte aus dem Zylinder herauskommt.
Und die rot-rote Karte wiederum präsentiert dem Laborassistenten in
100% aller Fälle eine rote Seite, also setzt Zeile 13 den Multiplikator auf
den Wert 1.
Abschließend gibt der Aufruf der Methode prob("RR")
die normalisierte
Wahrscheinlichkeit für die Hypothese "Rot-rote Karte" zurück, die
2/3 beträgt.
01 #!/usr/local/bin/perl -w 02 use strict; 03 use Distrib; 04 05 my $distrib = Distrib->new(); 06 07 $distrib->set( "BB", 0.33 ); 08 $distrib->set( "BR", 0.33 ); 09 $distrib->set( "RR", 0.33 ); 10 11 $distrib->mult( "BB", 0 ); 12 $distrib->mult( "BR", 0.5 ); 13 $distrib->mult( "RR", 1 ); 14 15 $distrib->normalize(); 16 17 print $distrib->prob( "RR" ), "\n";
Listing 3 nutzt das CPAN-Modul Moose, um sich die Definition des Distrib-Konstruktors in Perl zu ersparen, allerdings muss der Code statt dessen mit "has" die Klassenattribute deklarieren und initialisieren. Im vorliegenden Fall nur ein netter Gimmick, bei weiteren Attributen wäre der Code deutlich schlanker als bei manueller Klassenfahrt.
01 use Moose; 02 03 has 'values' => (is => 'rw', isa => 'HashRef', 04 default => sub { {} } ); 05 06 sub set { 07 my( $self, $hypo, $prob ) = @_; 08 09 $self->values()->{ $hypo } = $prob; 10 } 11 12 sub mult { 13 my( $self, $hypo, $prob ) = @_; 14 15 $self->values()->{ $hypo } *= $prob; 16 } 17 18 sub normalize { 19 my( $self ) = @_; 20 21 my $values = $self->values(); 22 my $sum = 0; 23 24 for my $hypo ( keys %$values ) { 25 $sum += $values->{ $hypo }; 26 } 27 for my $hypo ( keys %$values ) { 28 $values->{ $hypo } /= $sum; 29 } 30 } 31 32 sub prob { 33 my( $self, $hypo ) = @_; 34 35 return $self->values()->{ $hypo }; 36 } 37 38 1;
Die Methode set()
nimmt den Namen einer Hypothese (zum Beispiel "BB")
und deren A-priori-Wahrscheinlichkeit entgegen, und speichert sie in dem
Objekt-internen Array values
. Zum Multiplizieren eines Hypothesenwerts
mit einem konstanten Wert nimmt mult()
beide entgegen, findet den
bislang gespeicherten Wert in values
und multipliziert ihn mit dem
hereingereichten Wert für $prob
. Die Methode normalize()
iteriert
über alle bislang eingespeisten Werte im Hash, summiert sie in $sum
auf und dividiert dann alle Werte durch diesen Wert. So beträgt die
neue Summe aller Wahrscheinlichkeitswerte zum Beispiel nach einer
Multiplikation wieder 1, jeder Wert kann also als Wahrscheinlichkeit
zwischen 0 und 1 interpretiert werden. Am Ende findet prob()
den
Wert für eine gesuchte Hypothese, indem es den Werte-Hash mittels der
von Moose automatisch erzeugten Methode values()
als Referenz hervorholt
und unter dem Schlüssel der gesuchten Hypothese den dort abgelegten Wert
extrahiert.
Führt man mehrere solcher Tests für verschiedenartigste Probleme durch,
zeigt sich ein wiederkehrendes Muster: Nach dem Aufsetzen der Hypothesen
werden die Wahrscheinlichkeiten zusätzliche Daten jeweils mit allen
definierten Hypothesen multipliziert. Es bietet sich also wie in [2]
vorgemacht an, eine von Distrib
abgeleitete Musterklasse HypoTest
wie in Listing 2 zu definieren, die mittels einer Methode update()
die Werte für alle Hypothesen auffrischt. Weiter verlässt sie sich darauf,
dass wiederum von ihr abgeleitete Klassen (wie z.B. CardHypoTest
in Listing 4) die abstrakte Methode likelihood()
überladen und
basierend auf der übergebenen Hypothese und den zusätzlich verfügbaren
Daten die Wahrscheinlichkeit P(D|H) zurückgeben. Das Framework HypoTest
ruft dann diese Methode wiederholt auf, um die Werte für einzelne Hypothesen
einzuholen, bevor es diese in der Verteilung Distrib
setzt.
Das Framework bietet weiterhin eine Methode print()
, mit der es die
Werte aller aktualisierten Wahrscheinlichkeiten für die jeweilige Hypothese
ausgibt.
01 use Moose; 02 use Distrib; 03 04 has 'distrib' => ( 05 is => 'rw', 06 isa => 'Distrib', 07 default => sub { Distrib->new() } ); 08 09 sub hypo_add { 10 my( $self, $hypo ) = @_; 11 12 $self->distrib()->set( $hypo, 1 ); 13 } 14 15 sub update { 16 my( $self, $data ) = @_; 17 18 my $hypos = $self->distrib()->values(); 19 20 for my $hypo ( keys %$hypos ) { 21 $self->distrib()->mult( $hypo, 22 $self->likelihood( $data, $hypo ) 23 ); 24 } 25 26 $self->distrib()->normalize(); 27 } 28 29 sub print { 30 my( $self ) = @_; 31 32 my $values = $self->distrib()->values(); 33 34 for my $hypo ( keys %$values ) { 35 print "$hypo: $values->{ $hypo }\n"; 36 } 37 } 38 39 sub likelihood { 40 die "Subclass needs to override this"; 41 } 42 43 1;
So nimmt likelihood()
in Listing 4 in $data vom Hauptprogramm den
Buchstaben "R" entgegen (um eine gezogene Karte mit roter Vorderseite als
Zusatzbedingung anzumelden) und rechnet dann basierend auf der ebenfalls
hereingereichten Hypothese ("RR", "RB", "BB") aus, wie wahrscheinlich es
ist, dass der Proband auf eine rote Vorderseite blickt: Eine glatte
1 (also 100%) für "RR", 0.5 für "RB" und 0 für "BB".
Die Funktion setzt dazu einen regulären Ausdruck ab, der die Anzahl der
"R"s in der Hypothese zählt, und dividiert das Ergebnis anschließend durch
zwei als Fließkommawert, damit Perl keine Integerdivision vornimmt
und den Divisionsrest unterschlägt.
Am Ende gibt hypotest
die Wahrscheinlichkeiten für alle Hypothesen
in der Verteilung mittels der Methode print()
aus dem Modul HypoTest
an und berichtet zurecht, dass die rot-rote Karte
in 2/3 aller Fälle erscheint und damit der Proband, der auf eine
rote Vorderseite blickt, nach dem Umdrehen der Karte auch eine rote
Rückseite sieht:
$ ./hypotest RB: 0.333333333333333 RR: 0.666666666666667 BB: 0
01 #!/usr/local/bin/perl -w 02 use strict; 03 package CardHypoTest; 04 use base qw( HypoTest ); 05 06 sub likelihood { 07 my( $self, $data, $hypo ) = @_; 08 09 # count the number or Rs in the string 10 my @matches = ( $hypo =~ /($data)/g ); 11 return @matches / 2.0; 12 } 13 14 package main; 15 my $suite = CardHypoTest->new(); 16 17 $suite->hypo_add( "RR" ); 18 $suite->hypo_add( "BB" ); 19 $suite->hypo_add( "RB" ); 20 21 $suite->update( "R" ); 22 $suite->print();
Mit einfachen Skripts lassen sich so mathematisch fundierte Folgerungen aufstellen, und auch wer schon beim Anblick einer Formel in Ohnmacht fällt, kommt so auf seine Kosten und lernt, wie weitere Beispiele in [2] oder [3] (zwar in Python und nicht in Perl) zeigen, dass statistische Alltagsprobleme oft erstaunliche Lösungen haben, die auf den ersten Blick gar nicht intuitiv zu erfassen sind.
Listings zu diesem Artikel: ftp://www.linux-magazin.de/pub/listings/magazin/2014/07/Perl
"Think Bayes", Allen B. Downey, O'Reilly, 2013
"Probabilistic-Programming-and-Bayesian-Methods-for-Hackers", Cameron Davidson-Pilon, https://github.com/CamDavidsonPilon/Probabilistic-Programming-and-Bayesian-Methods-for-Hackers
"Ziegenproblem" - Statistische Probleme -- geknackt mit Perl, Michael Schilli, Linux-Magazin 09/1999
Satz von Bayes, http://de.wikipedia.org/wiki/Satz_von_Bayes
Ziegenproblem, http://de.wikipedia.org/wiki/Ziegenproblem