Echsen-Beispiel

In einem biologischen Experiment werden Echsen markiert und später nochmals eingefangen. Dabei ist erstmal von Interesse, ob eine markierte Echse wieder eingefangen werden kann.

Das Einfangen ist ein Bernoulli-Experiment, ein Erfolg ist dabei, dass die markierte Echse wieder eingefangen werden kann. Bezeichne \(p_t\) die Einfangwahrscheinlichkeit nach \(t\) Zeiteinheiten.

Für die uns unbekannte Einfangwahrscheinlichkeit \(p_t\) können die Biologen a priori eine gute Einschätzung geben. Sie nennen folgende Einschätzungen:

Zeitpunkt	2	3	4	5
Mittelwert	0.3	0.4	0.5	0.2
95% Kredibilitätsintervall	[0.1,0.5]	[0.2,0.6]	[0.3,0.7]	[0.05,0.4]

Gewinnung Subjektiver Prioris

Subjektive Prioris basieren auf Vorwissen eines (oder mehrerer) Experten über die Verteilung des Parameters. Das Vorwissen liegt dabei nicht unbedingt in Form einer Dichte vor.

Für diskrete Parameter ist Angabe der einzelnen Wahrscheinlichkeiten relativ einfach (unter der Vorgabe, dass sich die Wahrscheinlichkeiten aus eins addieren).

Für stetige Parameter bestehen folgende Möglichkeiten:

• Diskretisierung des Parameters. Wie beim Histogramm kann also die Wahrscheinlichkeit für mehrere Intervalle angegeben werden.
• Parametrische subjektive Prioris.
• Maximum Entropy Priori.

Als ersten Schritt geben wir uns die Form der Verteilung vor, zum Beispiel im einfachsten Fall die konjugierte Priori.

Ist die Form bekannt, können wir aus den Vorinformationen die Parameter bestimmen, z.B. aus dem Mittelwert oder die Varianz. Günstig ist es dabei, robuste Statistiken wie Median oder Quantil zu verwenden.

Schliesslich, soweit möglich, ist zu überprüfen, ob Vorinformationen und gewählte Form zusammenpassen. Wählt man zum Beispiel die Normalverteilung als Form, die Vorinformation weisen aber auf eine schiefe Verteilung hin, muss man die Verteilungsannahme korrigieren.

Im Echsenbeispiel

Für den Zeitpunkt \(t=2\) war die Vorinformation der Biologen:

• Mittelwert 0.3
• 95%-Intervall: [0.1,0.5]

Wir benutzen die Beta-Verteilung als konjugierte Verteilung. Der Erwartungswert der Beta-Verteilung ist \(\frac{a}{a+b}\). Es ergibt sich:

\[ \frac{a}{a+b} = 0.3 \Leftrightarrow b = \frac{7}{3} a \]

Dann müssen wir \(a\) so suchen, dass \(P(0.1<X<0.5)=0.95\), gilt, wobei \(X\sim \operatorname{Beta}\left(a,\tfrac{7}{3}a\right)\)ist. Zum Beispiel kann man in R die Funktion optimize() benutzen:

fehler<-function(a)
  {
  wahrscheinlichkeit <- pbeta(0.5,a,7*a/3)-pbeta(0.1,a,7*a/3)
  return((wahrscheinlichkeit-0.95)^2)
  }
a<-optimize(fehler,c(.01,100))
print(a$minimum)

## [1] 5.45238

• Wir definieren uns eine Funktion fehler, die die quadrierte Differenz zwischen 0.95 und der Wahrscheinlichkeit für das Intervall [0.1, 0.5] einer Beta\((a,\frac{7}{3}a)\)-Verteilung für einen gegebenen Wert a liefert
• pbeta(q,a,b) liefert den Wert der Verteilungsfunktion der Beta\((a,b)\)-Verteilung an der Stelle \(q\)
• optimize minimiert die Funktion fehler bezüglich a; c(.01,100) gibt dabei den Bereich an, in dem optimize sucht.

In diesem Beispiel ergibt sich als Vorinformation der Experten (gerundet):

\[ p_2 \sim \text{Beta}(5.45,12.72) \]

Weitere Zeitpunkte

Zeitpunkt	2	3	4	5	6
Mittelwert	0.3	0.4	0.5	0.2	0.2
95% Kredibilitätsintervall	[0.1,0.5]	[0.2,0.6]	[0.3,0.7]	[0.05,0.4]	[0.05,0.4]

• \(p_3 \sim\)Beta(\(8.58, 12.87)\)
• \(p_4 \sim\)Beta(\(11.26, 11.26)\)
• \(p_5 \sim\)Beta(\(3.50, 14.00)\)
• \(p_6 \sim\)Beta(\(3.50, 14.00)\)

Entropie

Die Shannon-Entropie (nach Claude Shannon, Begründer der Informationstheorie) ist für diskrete \(\theta\) definiert als

\[ H(\theta) = -\sum_\theta p(\theta)\log(p(\theta)) \]

Für stetige \(\theta\) gilt

\[ H(\theta) \propto -\int_\theta p(\theta)\log(p(\theta)) \]

Maximum-Entropie

Will man eine Priori mit möglichst wenig Information, so maximiert man die Entropie bezüglich der Dichte. Sind Momente \(E(g_k(X))\) bekannt, dann lässt sich eine Maximum-Entropie-Priori (MEP) in der Form

\[ p_{\text{MEP}}=c\exp\left(\sum_{k=1}^K\lambda_kg_k(x)\right) \]

finden, so dass die Momente erhalten bleiben.

Beispiel

Sei \(\theta\) ein unbekannter Parameter, für den wir eine Priori konstruieren wollen.

Der Experte sagt: A priori hat \(\theta\):

• \(E(\theta)=1\)
• Var\((\theta)=4\)

Aus dem Verschiebungssatz ergibt sich

\[ \text{Var}(\theta) = E(\theta^2)+\left(E(\theta)\right)^2 \\ \Rightarrow E(\theta^2) = \text{Var}(\theta) - \left(E(\theta)\right)^2 = 4-1 = 3 \]

Hier sind also die Momente von \(g_1(\theta)=\theta\) und \(g_2(\theta)=\theta^2\) bekannt.

Damit hat die Maximum-Entropy-Priori die Form

\[ p_{\text{MEP}} = c\cdot \exp(\lambda_1 g_1(\theta) + \lambda_2 g_2(\theta))\\ = c\cdot \exp(\lambda_1 \theta + \lambda_2 \theta^2) \]

Nun muss man \(\lambda_1\) und \(\lambda_2\) sowie die Normalisierungskonstante \(c\) so bestimmen, das es sich um eine Dichte handelt.

Falls \(\theta\) auf ganz \(\mathbb{R}\) definiert ist, kann man zeigen, dass es sich bei dieser Dichte um eine Normalverteilungsdichte handelt (dazu später mehr im Abschnitt Normalverteilungsmodell):

\[ p_{\text{MEP}} = \frac{1}{\sqrt{2\cdot\pi\cdot 4}} \exp\left(\frac{1}{2\cdot 4}(\theta-1)^2\right) \]

Die MEP ist dann also eine N\((1,4)\)-Verteilung.

Ist \(\theta\) beschränkt, kommt man hier zu einer trunkierten Normalverteilung.

In allgemeineren Fällen müssen die \(\lambda_i\) numerisch bestimmt werden. Insbesondere dann, wenn mehrere Priori-Informationen vorliegen.