Bestimmung von Parametern
Die Maximum Likelihood Methode
Gegeben sei eine zufällige Veränderliche $x$ mit der Dichtefunktion $p(a;x)$. $a$ sei ein zunächst unbekannter Parameter der Dichtefunktion. Die Veränderliche $x$ sei das Ergebnis eines Versuches oder einer Simulation. Wir führen den Versuch (Simulation) $N$-mal durch und erhalten $N$ Meßergebnisse $x_{1},x_{2},...,x_{N}$ der Veränderlichen $x$. Unsere Aufgabe ist jetzt, aus diesen $N$ Werten den Parameter $a$ zu bestimmen, und zwar mit Hilfe einer Schätzfunktion
\begin{displaymath}
a = f(x_{1},x_{2},..., x_{N}).
\end{displaymath} (1)

Im Hinblick auf spätere Anwendungen ist es günstiger, diese Funktion in impliziter Form zu schreiben:
\begin{displaymath}
F(a; x_{1},x_{2},...,x_{N}) = f(x_{1},x_{2},x_{3},....,x_{n}) - a = 0.
\end{displaymath} (2)

Die Meßergebnisse $x_{1},x_{2},...x_{N}$ nennen wir eine Stichprobe aus der Menge der Veränderlichen $x$. Der mit Hilfe der Gleichung (1) oder (2) bestimmte Parameter $a$ ist nun wieder eine zufällige Veränderliche, da er aus endlich vielen zufälligen Veränderlichen berechnet wurde. Der Parameter $a$ besitzt also eine Dichtefunktion, die wir im folgenden mit $q(a)$ bezeichnen werden. Die Veränderlichen $x_{1},x_{2},...,x_{N}$ sind unabhängige Veränderliche, da sie aus der mehrmaligen Wiederholung ein und desselben Versuches erhalten wurden. Wir können daher $\vec{x} = (x_{1},x_{2},...,x_{N})$ als Vektor auffassen, mit der Dichtefunktion
\begin{displaymath}
p(a;x_{1},x_{2},...,x_{N}) = \prod_{\nu=1}^{N} p(a; x_{\nu}).
\end{displaymath} (3)

Für eine Summe $y$ von unabhängigen Veränderlichen,

\begin{displaymath}
y = x_{1} + x_{2} + .... + x_{N}
\end{displaymath}

hatten wir bereits gezeigt, daß $y$ näherungsweise normalverteilt war, zumindest für große $N$. Dieses Ergebnis kann man verallgemeinern. Für nicht zu exotische Funktionen $F(a;x_{1},x_{2},...,x_{N})$ ist $a$ ebenfalls näherungsweise normalverteilt. Wir dürfen daher davon ausgehen, daß $q(a)$ durch einen Erwartungswert $<a>$ und Varianz $\sigma_{a}$ vollständig beschrieben werden kann. Das Problem ist dann, eine Funktion $F(a;x_{1},x_{2},...,x_{N})$ so zu bestimmen, daß die Varianz $\sigma_{a}$ der Dichtefunktion $q(a)$ möglichst klein wird. Ein solches Verfahren ergibt den kleinsten Meßfehler bei vorgegebener Stichprobe $x_{1},x_{2},...,x_{N}$. Die Antwort der Statistik auf diese Aufgabe ist das sogenannte Maximum- Likelihood Verfahren:
Satz: Die beste Schätzfunktion für die Bestimmung des Parameters $a$ ist die Ableitung des Logarithmus der Dichtefunktion nach dem Parameter $a$,

\begin{displaymath}
F(a; x_{1},x_{2},...,x_{N}) =
\frac{\partial ln p(a;x_{1},x_{2},...x_{N})}{\partial a}.
\end{displaymath} (4)

Die aus der Bestimmungsgleichung
\begin{displaymath}
F(a; x_{1},x_{2},...,x_{N}) = 0
\end{displaymath} (5)

berechnete zufällige Veränderliche $a$ ist näherungsweise normalverteilt mit der Varianz
\begin{displaymath}
\sigma_{a} = \left( -\langle \frac{\partial F}{\partial a} \...
...},x_{2},...,x_{N})}
{\partial^{2} a} \rangle \right)^{-1/2} .
\end{displaymath} (6)

Die Funktion
\begin{displaymath}
S(a) = ln p(a; x_{1},x_{2},...,x_{N})
\end{displaymath} (7)

nennt man die Likelihood Funktion.
Man kann also, und das ist für spätere Anwendungen von Bedeutung, das Maximum- Likelihood Verfahren als Extremalaufgabe schreiben:
\begin{displaymath}
S(a) =! \; \; Maximum.
\end{displaymath} (8)

Den aus einer einzigen Meßreihe $x_{1},x_{2},...,x_{N}$ gewonnenen Wert für $a$ nennen wir im folgenden den Schätzwert $\overline{a}$. Der Fehler für $\overline{a}$ kann näherungsweise aus
\begin{displaymath}
\Delta \overline{a} \approx \left( -\frac{\partial F(a; x_{1},x_{2},...,x_{N})}
{\partial a} \right)_{a=\overline{a}}^{-1/2}
\end{displaymath} (9)

berechnet werden.

Beispiel. Eine einparametrige Wahrscheinlichkeitsverteilung ist die Poisson Verteilung

\begin{displaymath}
P(a; k) = e^{-a} \frac{a^{k}}{k!}.
\end{displaymath}

Es seien $N$ Messungen $k_{\nu}$ der zufälligen Veränderlichen $k$ vorgelegt. Die Wahrscheinlichkeitsverteilung der unabhängigen Veränderlichen $k_{1},k_{2},...,k_{N}$ ist

\begin{displaymath}
P(a; k_{1},k_{2},...,k_{N}) = \prod_{\nu=1}^{N} e^{-a}
\frac{a^{k_{\nu}}}{k_{\nu}!}.
\end{displaymath}

Die Funktion $F$ wird zu

\begin{displaymath}
F(a; k_{1},k_{2},...,k_{N}) = -N + \frac{1}{a} \sum_{\nu=1}^{N} k_{\nu}.
\end{displaymath}

Aus $F=0$ erhalten wir den besten Schätzwert für $a$ zu

\begin{displaymath}
\overline{a} = \frac{1}{N} \sum_{\nu=1}^{N} k_{\nu} ,
\end{displaymath}

also den einfachen arithmetischen Mittelwert. Der Fehler der Messung berechnet sich aus

\begin{displaymath}
\frac{\partial F}{\partial a} = - \frac{1}{a^{2}} \sum_{\nu=1}^{N} k_{\nu}
= - \frac{1}{a^{2}} \overline{a} N
\end{displaymath}

zu

\begin{displaymath}
\Delta \overline{a} = \pm \left(
-\frac{\partial F}{\partia...
..._{a=\overline{a}}^{-1/2} =
\pm \sqrt{\frac{\overline{a}}{N}}.
\end{displaymath}

Mehrparametrige Verteilungen
Wir gehen weiter zu den Verteilungen mit mehreren Parametern:

\begin{displaymath}
p(\vec{a}; x) = p(a_{1},a_{2},...,a_{s}; x).
\end{displaymath} (10)

In diesem Fall sind $s$ Funktionen zu suchen,
\begin{displaymath}
F_{\sigma}(\vec{a}; x_{1},x_{2},...,x_{N}) = 0, \; \; \; \sigma=1,2,...,s,
\end{displaymath} (11)

und zwar so, daß die Norm der Kovarianzmatrix der Dichtefunktion $q(a_{1},a_{2},...,a_{s})$ möglichst klein wird. Auch für diesen Fall läßt sich zeigen, daß die Funktionen
\begin{displaymath}
F_{\sigma}(\vec{a};x_{1},x_{2},...,x_{N})
= \frac{\partial ...
.....,x_{N})}{\partial a_{\sigma}},
\; \; \; \sigma = 1,2,...,s,
\end{displaymath} (12)

die besten Schätzwerte für den Parametervektor $\vec{a}$ ergeben. Wie bei den einparametrigen Verteilungen nennen wir
\begin{displaymath}
S(\vec{a}) = ln p(\vec{a}; x_{1},x_{2},...,x_{N})
\end{displaymath} (13)

die Likelohood Funktion. Im Gegensatz zu (2) handelt es sich bei der Formel (11) allerdings um eine System mit $s$ Gleichungen und $s$ Unbekannten. Wir bilden die Ableitungen
\begin{displaymath}
\frac{\partial F_{\sigma}}{\partial a_{\tau}}
= \frac{\part...
...{1},x_{2},...,x_{N})}
{\partial a_{\sigma} \partial a_{\tau}}.
\end{displaymath} (14)

In Verallgemeinerung zu (9) läßt sich zeigen, daß die Fehlermatrix für die Schätzwerte $\overline{a}_{1}, \overline{a}_{2},...,\overline{a}_{s}$ durch
\begin{displaymath}
(\Delta \overline{a}_{\sigma} \Delta \overline{a}_{\tau}) =
...
...
{\partial a_{\tau}} \right)_{\vec{a}=\overline{\vec{a}}}^{-1}
\end{displaymath} (15)

gegeben ist. Man beachte, daß es sich hierbei um eine Matrizengleichung handelt, insbesondere auf der rechten Seite der Gleichung um eine Matrizeninversion. Als Fehler des Schätzwertes $\overline{a}_{\sigma}$ eines einzelnen Parameters gibt man im allgemeinen die Wurzel aus dem Diagonalelement $\Delta \overline{a}_{\sigma} \Delta \overline{a}_{\sigma}$ an.

Beispiel. Als Beispiel diskutieren wir die zweiparametrige Normalverteilung

\begin{displaymath}
p(a,\sigma ; x) = \frac{1}{\sqrt{2\pi} \sigma} e^{-(x-a)^{2}/(2 \sigma^{2})}.
\end{displaymath}

Es seien $N$ unabhängige Messungen $x_{1},x_{2},...,x_{N}$ der zufälligen Veränderlichen $x$ vorgegeben. Die gemeinsame Dichtefunktion ist

\begin{displaymath}
p(a, \sigma ; x_{1},x_{2},...,x_{N}) = \frac{1}{(2\pi)^{N/2}...
...ma^{N}}
\prod_{\nu=1}^{N} e^{-(x_{\nu}-a)^{2}/(2 \sigma^{2})}.
\end{displaymath}

Das Gleichungssystem zur Bestimmung der Schätzwerte $\overline{a}$ und $\overline{\sigma}$ wird zu
$\displaystyle \frac{1}{\sigma^{2}} \sum_{\nu=1}^{N} (x_{\nu}-a)$ $\textstyle =$ $\displaystyle 0 ,$  
$\displaystyle -\frac{N}{\sigma} + \frac{1}{\sigma^{3}} \sum_{\nu=1}^{N} (x_{\nu}-a)^{2}$ $\textstyle =$ $\displaystyle 0.$  

Die Lösung ist
$\displaystyle \overline{a}$ $\textstyle =$ $\displaystyle \frac{1}{N} \sum_{\nu=1}^{N} x_{\nu},$  
$\displaystyle \overline{\sigma}$ $\textstyle =$ $\displaystyle \sqrt{\frac{1}{N}
\sum_{\nu=1}^{N} (x_{\nu}-\overline{a})^{2}}.$  

Zur Bestimmung der Fehlermatrix berechnen wir
$\displaystyle \frac{\partial F_{a}}{\partial a}$ $\textstyle =$ $\displaystyle - \frac{N}{\sigma^{2}}$  
$\displaystyle \frac{\partial F_{a}}{\partial \sigma}$ $\textstyle =$ $\displaystyle - \frac{2}{\sigma^{3}}
\sum_{\nu=1}^{N} (x_{\nu}-a)$  
$\displaystyle \frac{\partial F_{\sigma}}{\partial a}$ $\textstyle =$ $\displaystyle - \frac{2}{\sigma^{3}}
\sum_{\nu=1}^{N} (x_{\nu}-a)$  
$\displaystyle \frac{\partial F_{\sigma}}{\partial \sigma}$ $\textstyle =$ $\displaystyle \frac{N}{\sigma^{2}}
- \frac{3}{\sigma^{4}} \sum_{\nu=1}^{N} (x_{\nu} -a)^{2}.$  

Wir bilden die Fehlermatrix an der Stelle $a=\overline{a}$ und $\sigma = \overline{\sigma}$ und erhalten

\begin{displaymath}
\left( \begin{array}{cc} \partial F_{a}/\partial a &
\parti...
...2} & 0 \\
0 & -2N/\overline{\sigma}^{2} \end{array} \right) .
\end{displaymath}

Die hierzu inverse Matrix ist

\begin{displaymath}
\left( \begin{array}{cc} \partial F_{a}/\partial a &
\parti...
...2}/N & 0 \\
0 & -\overline{\sigma}^{2}/2N \end{array} \right)
\end{displaymath}

Daraus ergeben sich die Fehler der Schätzwerte zu
$\displaystyle \Delta \overline{a}$ $\textstyle =$ $\displaystyle \pm \frac{\overline{\sigma}}{\sqrt{N}},$  
$\displaystyle \Delta \overline{\sigma}$ $\textstyle =$ $\displaystyle \pm \frac{\overline{\sigma}}{\sqrt{2N}}.$  

Mehrdimensionale Dichtefunktionen
Die Parameter Bestimmung bei mehrdimensionalen zufälligen Vektoren ergibt keine neue Situation. Der zufällige Vektor $\vec{x} = (x_{1},x_{2},...,x_{n})$ werde $N$-mal unabhängig voneinander gemessen, die Meßwerte seien

$\displaystyle \vec{x}_{1}$ $\textstyle =$ $\displaystyle (x_{11},x_{21},.....,x_{n1}),$  
$\displaystyle \vec{x}_{2}$ $\textstyle =$ $\displaystyle (x_{12},x_{22},.....,x_{n2}),$  
  $\textstyle .$    
  $\textstyle .$    
  $\textstyle .$    
$\displaystyle \vec{x}_{N}$ $\textstyle =$ $\displaystyle (x_{1N},x_{2N},.....,x_{nN}).$  

Die Kompenenten der Vektoren $\vec{x}_{\nu}$ müssen jetzt natürlich nicht unabhängig voneinander sein. Wir erhalten daher lediglich eine Zerlegung der gesamten Dichtefunktion in der Form
\begin{displaymath}
p(\vec{a}; \vec{x}_{1},\vec{x}_{2}, ..., \vec{x}_{N}) =
p(\v...
...}; \vec{x}_{N}) = \prod_{\nu=1}^{N} p(\vec{a}; \vec{x}_{\nu}).
\end{displaymath} (16)

Beispiel. Als Beispiel betrachten wir die mehrdimensionale Normalverteilung

\begin{displaymath}
p(A,\vec{a}; \vec{x}) = \frac{\sqrt{det A}}{(2\pi)^{-n/2}}
e...
... \mu =1}^{n} A_{\nu \mu} (x_{\nu}-a_{\nu})
(x_{\mu}-a_{\mu})}.
\end{displaymath} (17)

Wir nehmen an, daß die Elemente der Matrix $A$ sowie der Erwartungswerte $\vec{a}$ durch $s$ andere Parameter $b_{1},b_{2},...,b_{s}$ ausgedrückt werden können:
$\displaystyle A_{\nu \mu}$ $\textstyle =$ $\displaystyle A_{\nu \mu}(b_{1},b_{2},...,b_{s}), \; \; \; \nu,\mu=1,2,...,n,$ (18)
$\displaystyle a_{\nu}$ $\textstyle =$ $\displaystyle a_{\nu}(b_{1},b_{2},...,b_{s}), \; \; \; \nu=1,2,...,n.$ (19)

Dieses schließt offensichtlich den Fall ein, daß die Parameter $b_{\sigma}$ identisch mit den Parametern $A_{\nu \mu}$ und $a_{\nu}$ sind, und damit $s= n^{2} + n$.

Seien $N$ unabhängige Messungen des zufälligen Vektors $\vec{x}$ vorgelegt. Die gemeinsame Dichtefunktion ist

\begin{displaymath}
p(\vec{b}; \vec{x}_{1}, \vec{x}_{2},..., \vec{x}_{N}) =
C \p...
...=1}^{n} A_{\nu \mu}
(x_{\nu i}-a_{\nu})(x_{\mu i} - a_{\mu})}.
\end{displaymath} (20)

mit
\begin{displaymath}
C = \frac{(det A)^{N/2}}{(2\pi)^{Nn/2}}.
\end{displaymath} (21)

Zur Vereinfachung schreiben wir im folgenden den Ausdruck im Exponenten in Matrizen- und Vektorform (siehe Kap.3.2):
\begin{displaymath}
p(\vec{b};\vec{x}_{1},\vec{x}_{2},...,\vec{x}_{N}) =
C \prod...
...\frac{1}{2} (\vec{x}_{i}-\vec{a}) A
(\vec{x}_{i} - \vec{a})}.
\end{displaymath} (22)

Die Ableitungen vom Logarithmus dieser Dichtefunktion sind wie folgt:

\begin{displaymath}
\frac{\partial ln p}{\partial b_{\sigma}} =
\frac{\partial ...
..._{i}-\vec{a}) A
\frac{\partial \vec{a}}{\partial b_{\sigma}}.
\end{displaymath}

Wegen

\begin{displaymath}
0 = \int d\vec{x}_{1} d\vec{x}_{2} \cdot \cdot \cdot d\vec{x...
...ac{\partial A_{\nu\mu}}{\partial b_{\sigma}} (A^{-1})_{\nu\mu}
\end{displaymath}

folgt

\begin{displaymath}
\frac{\partial ln C}{\partial b_{\sigma}} = \frac{N}{2}
\sum...
...c{\partial A_{\nu\mu}}{\partial b_{\sigma}}
(A^{-1})_{\nu\mu}.
\end{displaymath}

Bei dieser Umformung haben wir benutzt, daß
$\displaystyle <x_{\nu i}-a_{\nu}>$ $\textstyle =$ $\displaystyle 0 ,$ (23)
$\displaystyle <(x_{\nu i}-a_{\nu})(x_{\mu i}-a_{\mu})>$ $\textstyle =$ $\displaystyle (A^{-1})_{\nu\mu},$ (24)
$\displaystyle \int d\vec{x}_{1} d\vec{x}_{2} \cdot \cdot \cdot d\vec{x}_{N}
p(\vec{b}; \vec{x}_{1}, \vec{x}_{2},..., \vec{x}_{N})$ $\textstyle =$ $\displaystyle 1.$ (25)

Zusammengefaßt erhalten wir
$\displaystyle F_{\sigma} = \frac{\partial ln p}{\partial b_{\sigma}} =$ $\textstyle +$ $\displaystyle \frac{N}{2}
\sum_{\nu\mu=1}^{n} \frac{\partial A_{\nu\mu}}{\partial b_{\sigma}}
(A^{-1})_{\nu\mu}$ (26)
  $\textstyle +$ $\displaystyle \sum_{i=1}^{N} (\vec{x}_{i}-\vec{a}) A
\frac{\partial \vec{a}}{\p...
...x}_{i} -\vec{a}) \frac{\partial A}{\partial b_{\sigma}}
(\vec{x}_{i}-\vec{a}) .$  

Die besten Parameterwerte nach dem Maximum Likelihood Verfahren findet man dann aus dem Gleichungssytem
\begin{displaymath}
F_{\sigma}(\vec{a}; \vec{x}_{1}, \vec{x}_{2}, ..., \vec{x}_{N}) = 0, \; \; \;
\sigma = 1,2,...,s.
\end{displaymath} (27)

Dieses ist ein hochgradig nichtlineares Gleichungssystem und kann im allgemeinen nur mit numerischen Methoden gelöst werden. Hierzu werden wir wir am Ende dieses Kapitels ein Verfahren angeben. Die zweiten Ableitungen sind einfach und ergeben
$\displaystyle \frac{\partial F_{\sigma}}{\partial b_{\tau}} =$ $\textstyle +$ $\displaystyle \frac{\partial^{2} ln p}{\partial b_{\sigma} \partial b_{\tau}}$ (28)
$\displaystyle =$ $\textstyle +$ $\displaystyle \frac{N}{2} \sum_{\nu\mu=1}^{n} \left[
\frac{\partial^{2} A_{\nu\...
...partial b_{\sigma}}
\frac{\partial(A^{-1})_{\nu\mu}}{\partial b_{\tau}} \right]$  
  $\textstyle +$ $\displaystyle \sum_{i=1}^{N} \left[ -\frac{\partial \vec{a}}{\partial b_{\tau}}...
...}) A
\frac{\partial^{2} \vec{a}}{\partial b_{\sigma} \partial b_{\tau}} \right]$  
  $\textstyle +$ $\displaystyle \sum_{i=1}^{N} \left[ (\vec{x}_{i}-\vec{a})
\frac{\partial A}{\pa...
...l^{2} A}{\partial b_{\sigma} \partial b_{\tau}}
(\vec{x}_{i}-\vec{a}) \right] .$  

Wir bilden die Erwartungswerte der zweiten Ableitungen. Mit Hilfe der Formel (4.xx) erhalten wir
\begin{displaymath}
-\langle \frac{\partial F_{\sigma}}{\partial b_{\tau}} \rang...
...l b_{\tau}} A
\frac{\partial \vec{a}}{\partial b_{\sigma}} .
\end{displaymath} (29)

Wegen der allgemeinen Matrizengleichung (siehe (3.xx))

\begin{displaymath}
\frac{\partial A}{\partial b} = - A \frac{\partial A^{-1}}{\partial b} A
\end{displaymath}

und mit Einführung der Spur kann man diesen Ausdruck auch in der kompakten Form
\begin{displaymath}
\left(-\langle \frac{\partial F_{\sigma}}{\partial b_{\tau}}...
...ial b_{\tau}} A
\frac{\partial \vec{a}}{\partial b_{\sigma}}
\end{displaymath} (30)

schreiben. Die inverse Matrix hierzu ist dann die Fehlermatrix. Wir haben dieses Beispiel in voller Länge durchgerechnet, da es uns in den Anwendungen noch häufiger beschäftigen wird.

Die Näherung der kleinsten Quadrate
Das soeben durchgerechnete Maximum Likelihood Verfahren für die $n$-dimensionale Normalverteilung kann stark vereinfacht werden, sobald die Kovarianzmatrix bekannt ist, und somit nicht von den Parametern $\vec{b}$ abhängt. In diesem Fall ist

\begin{displaymath}
\frac{\partial A_{\nu\mu}}{\partial b_{\sigma}} = 0, \; \; \;
\sigma = 1,2,...,s.
\end{displaymath} (31)

Die Gleichungen (4.xx) und (4.xx) vereinfachen sich zu
$\displaystyle F_{\sigma}$ $\textstyle =$ $\displaystyle \sum_{i=1}^{N} (\vec{x}_{i}-\vec{a}) A
\frac{\partial \vec{a}}{\partial b_{\sigma}},$ (32)
$\displaystyle \left( -
\langle \frac{\partial F_{\sigma}}{\partial b_{\tau}} \rangle \right)$ $\textstyle =$ $\displaystyle N \frac{\partial \vec{a}}{\partial b_{\tau}} A
\frac{\partial \vec{a}}{\partial b_{\sigma}} .$ (33)

Die Schätzfunktionen $F_{\sigma}$ hätte man in diesem Fall auch aus der Extremalbedingung
\begin{displaymath}
\chi^{2} = \sum_{i=1}^{N} (\vec{x}_{i} - \vec{a}) A (\vec{x}_{i}-\vec{a})
=! \; \; Minimum
\end{displaymath} (34)

schließen können. Falls $A_{\nu\mu} = \delta_{\nu\mu}/\sigma_{\nu}^{2}$ eine Diagonalmatrix ist, vereinfacht sich (4.xx) weiter zu
\begin{displaymath}
\chi^{2} = \sum_{i=1}^{N} \sum_{\nu=1}^{n} \frac{1}{\sigma_{\nu}^{2}}
(x_{i\nu} - a_{i\nu})^{2} =! \; \; Minimum.
\end{displaymath} (35)

In dieser letzten Form nennt man das Verfahren auch die Methode der kleinsten Quadrate. Dieser Methode wollen wir uns im folgenden näher zuwenden.

Die Methode der kleinsten Quadrate
Zur Einführung in die Methode der kleinsten Quadrate diskutieren wir zunächst den Fall einer diskreten Veränderlichen $k$ mit der Wahrscheinlichkeitsverteilung $P(\vec{a}; k)$. $\vec{a}$ sei ein Vektor von $s$ unbekannten Parametern,

\begin{displaymath}
\vec{a} = (a_{1}, a_{2}, ..., a_{s}).
\end{displaymath}

Es seien $M$ unabhängige Messungen der Veränderlichen $k$ vorgelegt, und zwar $k_{1}, k_{2}, ..., k_{M}$. Für hinreichend große $M$ sind natürlich nicht alle Messungen $k_{\nu}$ voneinander verschieden. Wir nehmen an, daß ein bestimmter Wert $k_{\mu}$ mit einer Häufigkeit $l_{\mu}$ in der Meßreihe auftritt und setzen die Anzahl der verschiedenen gemessenen Werte der Veränderlichen $k$ gleich $N$. Dann ist

\begin{displaymath}
\sum_{\nu=1}^{N} l_{\nu} = M.
\end{displaymath}

Die relative Häufigkeit $h_{\nu} = l_{\nu}/M$ sollte dann, wie wir aus Kap.2.4 her wissen, für große $M$ gegen die Wahrscheinlichkeit $P(\vec{a}; k_{\nu})$ konvergieren:
\begin{displaymath}
\lim_{M \to \infty} h_{\nu} = P(\vec{a}; k_{\nu}).
\end{displaymath} (36)

Daraus ergibt sich offensichtlich die folgende Bedingung für die Bestimmung der Parameter $\vec{a}$:
Satz: Bei der Methode der kleinsten Quadrate erhält man die besten Parameterwerte aus der Extremalbedingung

\begin{displaymath}
\chi^{2} = \sum_{\nu=1}^{N} w_{\nu} (P(\vec{a}; k_{\nu})-h_{\nu})^{2}
=! \; \; Minimum.
\end{displaymath} (37)

Das Gewicht $w_{\nu}$ ist hierbei durch das inverse Quadrat des Fehlers der relativen Häufigkeit gegeben:
\begin{displaymath}
w_{\nu} = \frac{1}{(\Delta h_{\nu})^{2}}.
\end{displaymath} (38)

Für $\Delta h_{\nu}$ kann man entweder die Wurzel aus $h_{\nu}$ oder die Wurzel aus $P(\vec{a}; k_{\nu})$ wählen, d.h.
\begin{displaymath}
\Delta h_{\nu} = \sqrt{h_{\nu}},
\end{displaymath} (39)

oder
\begin{displaymath}
\Delta h_{\nu} = \sqrt{P(\vec{a}; h_{\nu})}.
\end{displaymath} (40)

Im Grenzfall $M \to \infty$ ergeben beide Fehler dasselbe Ergebnis. Im letzteren Fall (40) tritt natürlich das Problem auf, daß man die Parameterwerte $\vec{a}$ bereits kennen muß, um $\Delta h_{\nu}$ zu berechnen. Man behilft sich hierbei, indem man möglichst gute, bereits bekannte Werte für die Parameter in die Formel einsetzt. Im übrigen gibt es seit Einführung der Methode der kleinsten Quadrate einen ewigen Streit der Statistiker, welche der beiden Gewichtsfaktoren (39) oder (40) die bessere Parameterbestimmung ergibt. Einer Diskussion der Argumente wollen wir uns in diesem Buch nicht zuwenden.

Das Minimum für die $\chi^{2}$- Funktion wird an der Nullstelle der Ableitung angenommen,

\begin{displaymath}
\frac{\partial \chi^{2}}{\partial a_{\sigma}} = 2 \sum_{\nu=...
...;k_{\nu})-h_{\nu}) \frac{\partial P}{\partial a_{\sigma}} = 0.
\end{displaymath} (41)

Diese Bedingung ergibt die folgende Schätzfunktion für die Bestimmung der Parameter $\vec{a}$:
\begin{displaymath}
F_{\sigma} = \sum_{\nu=1}^{N} w_{\nu} (P(\vec{a}; k_{\nu})- h_{\nu})
\frac{\partial P}{\partial a_{\sigma}} = 0.
\end{displaymath} (42)

Die Fehler für die Schätzwerte $\overline{a}_{\sigma}, \sigma=1,2,...,s$ ergeben sich wiederum, wie wir es schon vom Maximum Likelihood Verfahren her kennen, aus den partiellen Ableitungen der Schätzfunktion $F_{\sigma}$ zu
\begin{displaymath}
(\Delta \overline{a}_{\sigma} \Delta \overline{a}_{\tau}) = ...
...{\partial a_{\tau}}
\right)_{\vec{a}=\overline{\vec{a}}}^{-1}.
\end{displaymath} (43)

Das Verfahren der kleinsten Quadrate kann ohne Schwierigkeit auf kontinuierliche Veränderliche $x$ mit der Dichtefunktion $p(\vec{a};x)$ übertragen werden. Hierzu führen wir eine Diskretisierung der Veränderlichen $x$ ein, d.h. wir bilden ein Histogramm und identifizieren
\begin{displaymath}
k_{\nu} \equiv \{ x; x \in [x_{\nu},x_{\nu +1}] \} .
\end{displaymath} (44)

Die Wahrscheinlichkeitsverteilung ist dann
\begin{displaymath}
P(\vec{a}; k_{\nu}) = \int_{x_{\nu}}^{x_{\nu +1}} dx p(\vec{a};x) .
\end{displaymath} (45)

In der allgemeinsten Form der Methode der kleinsten Quadrate sind $N$ unabhängige Messungen $y_{\nu}$ bei den Koordinaten $x_{\nu}$ vorgegeben. Die Fehler der Messungen seien mit $d y_{\nu}$ bezeichnet. Die Meßwerte $y_{\nu}$ sollen durch eine Funktion $y = f(\vec{a};x)$ mit unbekannten Parametern $a_{\sigma}, \sigma=1,2,...,s$ beschrieben werden. Hierbei braucht es sich nicht um zufällige Veränderliche zu handeln. Die $\chi^{2}$- Funktion wird definiert durch
\begin{displaymath}
\chi^{2} = \sum_{\nu=1}^{N} w_{\nu} (f(\vec{a};x_{\nu})-y_{\nu})^{2},
\end{displaymath} (46)

mit
\begin{displaymath}
w_{\nu} = \frac{1}{(dy_{\nu})^{2}}.
\end{displaymath} (47)

Die besten Schätzwerte für die Parameter $\vec{a}$ erhält man für das Minimum der $\chi^{2}$- Funktion,
\begin{displaymath}
\frac{\partial \chi^{2}}{\partial a_{\sigma}} = 0, \; \; \; \sigma=1,2,...s,
\end{displaymath} (48)

d.h. mit Hilfe der Schätzfunktionen
\begin{displaymath}
F_{\sigma} = \sum_{\nu=1}^{N} w_{\nu} (f(\vec{a}; x_{\nu})-y...
...artial f}{\partial a_{\sigma}} = 0, \; \; \; \sigma=1,2,...,s.
\end{displaymath} (49)

Die Fehler dieser Parameterbestimmung berechnet man wie bisher mit Hilfe der Matrizengleichung
\begin{displaymath}
(\Delta \overline{a}_{\sigma} \Delta \overline{a}_{\tau} ) =...
...\partial a_{\tau}}
\right)_{\vec{a}=\overline{\vec{a}}}^{-1}.
\end{displaymath} (50)

Beispiel. Im folgenden diskutieren wir eine typische Anwendung der Methode der kleinsten Quadrate. Es seien $N$ Messungen $y_{\nu}$ bei den Koordinaten $x_{\nu}$ vorgegeben. Es werde angenommen, daß die Meßwerte $y_{\nu}$ linear mit den Koordinaten $x_{\nu}$ zusammenhängen, d.h.

\begin{displaymath}
y = f(a,b; x) = ax + b.
\end{displaymath}

Die $\chi^{2}$- Funktion lautet

\begin{displaymath}
\chi^{2} = \sum_{\nu=1}^{N} w_{\nu} (a x_{\nu} + b - y_{\nu})^{2},
\end{displaymath}

und die besten Schätzwerte für $a$ und $b$ erhält man aus den Gleichungen
$\displaystyle F_{a}$ $\textstyle =$ $\displaystyle \sum_{\nu=1}^{N} w_{\nu} (a x_{\nu} + b - y_{\nu}) x_{\nu} = 0 ,$  
$\displaystyle F_{b}$ $\textstyle =$ $\displaystyle \sum_{\nu=1}^{N} w_{\nu} (a x_{\nu} + b - y_{\nu}) = 0 ,$  

oder, nach Auflösen der Summen,
$\displaystyle \overline{a} \sum_{\nu=1}^{N} w_{\nu} x_{\nu}^{2} + \overline{b}
\sum_{\nu=1}^{N} w_{\nu} x_{\nu}$ $\textstyle =$ $\displaystyle \sum_{\nu=1}^{N} w_{\nu} x_{\nu} y_{\nu}$  
$\displaystyle \overline{a} \sum_{\nu=1}^{N} w_{\nu} x_{\nu} + \overline{b}
\sum_{\nu=1}^{N} w_{\nu}$ $\textstyle =$ $\displaystyle \sum_{\nu=1} w_{\nu} y_{\nu}$  

Dieses Gleichungssystem läßt sich leicht nach $\overline{a}$ und $\overline{b}$ auflösen.

In Verallgemeinerung werde angenommen, daß die Meßwerte $y_{\nu}$ nicht linear mit den Koordinaten $x_{\nu}$ zusammenhängen, sondern daß diese Abhängigkeit durch ein Polynom $s$-ten Grades beschrieben werden kann,

\begin{displaymath}
y = a_{0} + a_{1} x + a_{2} x^{2} + ... + a_{s} x^{s} = \sum_{\mu=0}^{s}
a_{\mu} x^{\mu}.
\end{displaymath}

Wir erhalten dann

\begin{displaymath}
\chi^{2} = \sum_{\nu=1}^{N} w_{\nu} (\sum_{\mu=0}^{s}
a_{\mu} x_{\nu}^{\mu} - y_{\nu})^{2}
\end{displaymath}

und

\begin{displaymath}
F_{\sigma} = \sum_{\nu=1}^{N} w_{\nu} ( \sum_{\mu=0}^{s} a_{...
... - y_{\nu}) x_{\nu}^{\sigma} = 0, \; \; \; \sigma=0,1,2,...,s.
\end{displaymath}

Dieses ergibt ein lineares Gleichungssystem mit $s+1$ zu bestimmenden Schätzwerten $\overline{a}_{\mu}, \mu=0,1,2,...,s$:

\begin{displaymath}
\sum_{\mu=0}^{s} \overline{a}_{\mu} \sum_{\nu=1}^{N} w_{\nu}...
...{\nu} y_{\nu} x_{\nu}^{\sigma},
\; \; \; \sigma = 0,1,2,...,s.
\end{displaymath}

Die Poisson Methode
Eine etwas andere Methode wird häufig bei der Bestimmung von Parametern aus Messungen von Zählraten angewandt. Angenommen, wir messen Zählraten $z_{\nu}$ bei Koordinaten $k_{\nu}$, und erwarten aufgrund theoretischer Überlegungen Zählraten $Z(\vec{a}; k_{\nu})$. Dann ist die Wahrscheinlichkeit für eine Messung der Zählrate $z_{\nu}$ durch

\begin{displaymath}
P(\vec{a}; z_{\nu}) = \frac{Z(\vec{a};k_{\nu})^{z_{\nu}}}{z_{\nu}!}
e^{-Z(\vec{a};k_{\nu})}
\end{displaymath} (51)

gegeben. Bei $N$ unabhängigen Messungen $z_{\nu}, \nu=1,2,...,N$ bei den Koordinaten $k_{\nu}, \nu=1,2,...,N$ kann man die Gesamtwahrscheinlichkeit als Produkt schreiben,
\begin{displaymath}
P(\vec{a},k_{1},k_{2},...,k_{N}) = \prod_{\nu=1}^{N}
\frac{Z(\vec{a};k_{\nu})^{z_{\nu}}}{z_{\nu}!} e^{-Z(\vec{a};k_{\nu})}
\end{displaymath} (52)

und die besten Schätzfunktionen nach der Maximum Likelihood Methode sind
\begin{displaymath}
F_{\sigma} = \frac{\partial ln P}{\partial a_{\sigma}} =
\s...
...)
\frac{\partial Z(\vec{a};k_{\nu})}{\partial a_{\sigma}} = 0.
\end{displaymath} (53)

Auf dieses Verfahren werden wir in späteren Anwendungen noch zurückkommen.

Ein Rechnerprogramm zur Bestimmung von Parametern
Alle drei Methoden, die Maximum Likelihood Methode, die Methode der kleinsten Quadrate und die Poisson Methode können als Extremalaufgabe geschrieben werden:
1. Maximum Likelihood Methode:

\begin{displaymath}
L(\vec{a}) = \sum_{i=1}^{N} ln f(\vec{a}; x_{i}) =! \; \; \; Maximum.
\end{displaymath} (54)

2. Methode der kleinsten Quadrate:
\begin{displaymath}
\chi^{2}(\vec{a}) = \sum_{i=1}^{N} w_{i} (f(\vec{a};x_{i}) - y_{i})^{2} =!
\; \; \; Minimum.
\end{displaymath} (55)

3. Poisson- Methode:
\begin{displaymath}
P(\vec{a}) = \sum_{i=1}^{N} [ y_{i} ln f(\vec{a}; x_{i}) - f(\vec{a}; x_{i})]
=! \; \; \; Maximum.
\end{displaymath} (56)

Man beachte allerdings die verschiedene Bedeutung der Funktion $f(\vec{a};x)$ in diesen drei Ausdrücken. Die Bestimmungsgleichungen für die Parameter $a_{\sigma}, \; \sigma=1,2,...,s$ sowie die Ableitungen für die Berechnung der Fehler sind:
1. Maximum Likelihood Methode:

$\displaystyle F_{\sigma}$ $\textstyle =$ $\displaystyle \sum_{i=1}^{N} \frac{1}{f(\vec{a};x_{i})}
\frac{\partial f(\vec{a}; x_{i})}{\partial a_{\sigma}} = 0,$ (57)
$\displaystyle -\frac{\partial F_{\sigma}}{\partial a_{\tau}}$ $\textstyle =$ $\displaystyle \sum_{i=1}^{N}
\frac{1}{f^{2}(\vec{a};x_{i})} \left[ \frac{\parti...
...};x_{i})
\frac{\partial^{2} f}{\partial a_{\sigma} \partial a_{\tau}} \right] .$ (58)

2. Methode der kleinsten Quadrate:
$\displaystyle F_{\sigma}$ $\textstyle =$ $\displaystyle \sum_{i=1}^{N} w_{i} (y_{i}-f(\vec{a};x_{i}))
\frac{\partial f(\vec{a};x_{i})}{\partial a_{\sigma}} = 0,$ (59)
$\displaystyle -\frac{\partial F_{\sigma}}{\partial a_{\tau}}$ $\textstyle =$ $\displaystyle \sum_{i=1}^{N}
w_{i} \left[ \frac{\partial f}{\partial a_{\sigma}...
...;x_{i}))
\frac{\partial^{2} f}{\partial a_{\sigma} \partial a_{\tau}} \right] .$ (60)

3. Poisson- Methode:
$\displaystyle F_{\sigma}$ $\textstyle =$ $\displaystyle \sum_{i=1}^{N} (\frac{y_{i}}{f(\vec{a};x_{i})} -1)
\frac{\partial f(\vec{a};x_{i})}{\partial a_{\sigma}} = 0,$ (61)
$\displaystyle -\frac{\partial F_{\sigma}}{\partial a_{\tau}}$ $\textstyle =$ $\displaystyle \sum_{i=1}^{N}
\left( \frac{y_{i}}{f^{2}(\vec{a};x_{i})} \left[
\...
...right] - \frac{\partial^{2} f}{\partial a_{\sigma} \partial a_{\tau}} \right) .$ (62)

Zur Lösung dieser Gleichungen verwenden wir das einfache Newton Verfahren. Dieses beruht auf der folgenden Iteration, die wir zunächst für den Fall eines einzelnen Parameters erläutern. Beginnend mit einem Startwert $a^{(0)}$ erhält man einen besseren Wert für die Nullstelle der Funktion $F(a)$ gemäß
\begin{displaymath}
a^{(1)} = a^{(0)} - \left( \frac{\partial F}{\partial a}
\right)_{a=a^{(0)}}^{-1} F(a^{(0)}).
\end{displaymath} (63)

Die Iteration besteht aus der weiteren Anwendung dieser Rekursion,
\begin{displaymath}
a^{(i)} = a^{(i-1)} - \left( \frac{\partial F}{\partial a}
\right)_{a=a^{(i-1)}}^{-1} F(a^{(i-1)}),
\end{displaymath} (64)

bis entweder die Differenz $\vert a^{(i)}-a^{(i-1)}\vert$ oder der Betrag von $F(a^{(i)})$ einen bestimmten vorgegebenen Wert unterschreitet. Offensichtlich kann die Bedingung
\begin{displaymath}
\vert a^{(i)} - a^{(i-1)}\vert \ll \vert\Delta \overline{a}\...
...ial a} \right)_{a=a^{(i-1)}}^{-1}
\right\vert
\; F(a^{(i-1)})
\end{displaymath} (65)

zum Abbruch der Iteration herangezogen werden. Diese Iteration kann leicht auf Gleichungssysteme übertragen werden:
\begin{displaymath}
\vec{a}^{(i)} = \vec{a}^{(i-1)} - \left(
\frac{\partial F_{...
...ight)_{\vec{a}=\vec{a}^{(i-1)}}^{-1} \vec{F}(\vec{a}^{(i-1)}).
\end{displaymath} (66)

Die Fehlermatrix ist am Ende der Iteration direkt durch
\begin{displaymath}
(\Delta \overline{a}_{\sigma} \Delta \overline{a}_{\tau} ) \...
...}}{\partial a_{\tau}}
\right)_{\vec{a}=\vec{a}^{(i-1)}}^{-1}
\end{displaymath} (67)

gegeben, muß also nicht extra berechnet werden.

Ein Demonstrationsprogramm ist in dem folgenden Applet gegeben. Mit Hilfe eines 1. Generators erzeugen wir Daten mit einer bestimmten Verteilungsfunktion. Diese erste Verteilung kann noch mit einer 2. Verteilung gefaltet werden, sodaß man insgesamt eine grosse Menge verschiedener Verteilungsfunktionen erzeugen kann. Diese Gesamtverteilung kann dann mit dem Modell einer Normalverteilung, einer Exponentialverteilung oder einer Poissonverteilung gefitted werden. Als Beispiel zeigen wir in Abbildung 1 den Fit einer normalverteilten Simulation (linkes Fesnter) mit einer Normalverteilung (rechtes Fenster). Die gefitten Parameter (a1 = Mittelwert, a2 = Varianz) sind im Monitor- Fenster angezeigt.


Abbildung 1: Maximum Likelihood Fit einer normalverteilten Größe mit einer Normalverteilung.

In Abbildung 2 haben wir eine Bernoulli- Verteilung mit einer Exponentialverteilung gefaltet. Daher wird der Ausläufer zu hohen Werten der Veränderlichen nochmals stark erhöht. Diese Verteilung wurde ebenfalls mit einer Normalverteilung gefitted. Die Ergebnisse im Monitorfenster zeigen ebenfalls gute Ergebnisse, die simulierte Verteilung und die gefittete Verteilung sind jedoch total verschieden. Verteilung und Modell stimmen in diesem Fall nicht überein. Wie müssen daher noch ein Kriterium diskutieren, um die Übereinstimmung von experimentellen Daten und Modell zu prüfen.


Abbildung 2: Maximum Likelihood Fit einer verzerrten Bernoulliverteilung mit einem normalverteilten Modell.

Testen von Hypothesen
In vorigen Abschnitt haben wir diskutiert, wie die Parameter $\vec{a}$ einer Wahrscheinlichkeitsverteilung $P(\vec{a};\vec{k})$ bzw einer Dichtefunktion $p(\vec{a};\vec{x})$ aus $N$ Messungen der Veränderlichen $\vec{k}$ bzw $\vec{x}$ bestimmt werden können. Hierbei haben wir stillschweigend angenommen, daß die gemessenen Veränderlichen wirklich aus der Wahrscheinlichkeitsverteilung $P(\vec{a};\vec{k})$ bzw aus der Dichtefunktion $p(\vec{a};\vec{x})$ stammen. Diese Annahme ist in der Praxis jedoch nur selten gerechtfertigt. Die Veränderlichen $\vec{k}$ bzw $\vec{x}$ sind das Ergebnis eines Versuches mit unbekanntem Bedingungskomplex, die Wahrscheinlichkeitsverteilung bzw Dichtefunktion ist daher eine hypothetische Annahme. Es bleibt also noch zu prüfen, wie gut diese Hypothese mit den Messungen des Experiment übereinstimmt. Es ergibt sich daher folgende Situation. Wir haben $N$ Messungen der Veränderlichen $x$ (wir beschränken uns zunächst auf eindimensionale Veränderliche). Aus Überlegungen irgendwelcher Art stellen wir die Hypothese auf, daß die Dichtefunktion $p(\vec{a};x)$ den Bedingungskomplex unseres Versuches beschreibt und bestimmen mit einem der im vorigen Abschnitt diskutierten Verfahren die optimalen Werte der Parameter $\vec{a}$. Hierbei erhalten wir immer irdendwelche Parameterwerte, unabhängig davon, ob die Dichtefunktion $p(\vec{a};x)$ wirklich den Bedingungskomplex des Experiments beschreibt oder auch nicht. Dieses erkennt man im allgemeinen auch nicht an den Fehlern der Parameter. Die Fehler können beliebig klein sein, obwohl die Dichtefunktion völlig falsch war. Abschließend müssen wir daher testen, ob die Hypothese $p(\vec{a};x)$ für die vorgegebenen Messungen $x_{1},x_{2},...,x_{N}$ akzeptiert werden kann.

Hierzu gibt es sehr viele verschiedene Testverfahren, die zumindest zum Teil auf spezielle Probleme zugeschnitten sind. Wir hatten bereits die speziellen Verfahren des Korrelationstests und des Phasentests im Zusammenhang mit den Pseudozufallszahlen diskutiert. Allen Testverfahren gemeinsam ist, daß, ähnlich wie bei der Parameterbestimmung, eine skalare Testfunktion definiert wird,

\begin{displaymath}
f = f(x_{1},x_{2},...,x_{N}).
\end{displaymath} (68)

Die Größe $f$ ist wiederum eine zufällige Veränderliche mit einer Dichtefunktion $q(f)$. Diese Dichtefunktion kann durch die Dichtefunktion $p(\vec{a};x)$ ausgedrückt werden gemäß:
\begin{displaymath}
q(f) = \int_{f(x_{1},x_{2},...,x_{N})} dx_{1} dx_{2} \cdot \cdot \cdot
dx_{N} p(\vec{a};x_{1},x_{2},...,x_{N}),
\end{displaymath} (69)

mit

\begin{displaymath}
p(\vec{a};x_{1},x_{2},...,x_{N}) = \prod_{\nu=1}^{N} p(\vec{a};x_{\nu}).
\end{displaymath}

Ähnlich wie bei der Parameterbestimmung kann man auch hier nach einer ''besten'' Testfunktion fragen. Dieser Frage werden wir uns in diesem Buch jedoch nicht widmen, sondern sofort zu dem am häufigsten angewendeten Testverfahren übergehen, nämlich zum $\chi^{2}$- Test.

Die $\chi^{2}$- Verteilung
Wir sortieren die $N$ Messungen $x_{1},x_{2},...x_{N}$ in ein Histogramm mit $M$ Intervallen zwischen den Stützstellen $x^{(0)},x^{(1)},...,x^{(M)}$. Falls die hypothetisch angenommene Dichtefunktion $p(\vec{a};x)$ auf der gesamten reellen Achse definiert ist, setzen wir die unterste und oberste Intervallgrenze zu

\begin{displaymath}
x^{(0)} = -\infty, \; \; \; \; \; x^{(M)} = +\infty .
\end{displaymath}

Die Anzahl der Meßwerte in jedem dieser Intervalle sei $K_{1},K_{2},...,K_{N}$. Die Wahrscheinlichkeitsverteilung ist durch
\begin{displaymath}
P_{\mu} = \int_{x^{(\mu-1)}}^{x^{(\mu)}} dx p(\vec{a};x), \; \; \;
\mu = 1,2,...,M,
\end{displaymath} (70)

für kontinuierliche Veränderliche bzw durch
\begin{displaymath}
P_{\mu} = \sum_{x^{(\mu-1)} < k \leq x^{(\mu)}} P(\vec{a};k), \; \; \;
\mu = 1,2,...,M
\end{displaymath} (71)

für diskrete Veränderliche gegeben. Die hier geschilderte Situation entspricht nun genau der Definition der Multinominalverteilung. Daher ist die Wahrscheinlichkeit, genau $K_{\mu}$ $(\mu=1,2,...,M)$ Messungen in den Intervallen $(x^{(\mu-1)},x^{(\mu)}]$ zu erhalten, durch die Multinominalverteilung
\begin{displaymath}
P(K_{1},K_{2},...,K_{M-1}) = \frac{N!}{K_{1}! K_{2}! \cdot \...
...!}
P_{1}^{K_{1}} P_{2}^{k_{2}} \cdot \cdot \cdot P_{M}^{K_{M}}
\end{displaymath} (72)

mit
$\displaystyle P_{M}$ $\textstyle =$ $\displaystyle 1 - \sum_{\mu=1}^{M-1} P_{\mu},$ (73)
$\displaystyle K_{M}$ $\textstyle =$ $\displaystyle N - \sum_{\mu=1}^{M-1} K_{\mu}$ (74)

gegeben. Man beachte, daß diese Verteilung nur $M-1$ zufällige Veränderliche hat. Genau wie im Fall der einfachen Bernoulli- Verteilung (mit $M=2$) kann die Multinominalverteilung für große $N$ durch eine $(M-1)$- dimensionale Normalverteilung approximiert werden. Aus der erzeugenden Funktion der Multinominalverteilung,
\begin{displaymath}
M_{K_{1},K_{2},...,K_{M-1}}(v_{1},v_{2},...,v_{M-1}) =
[ 1+ \sum_{\mu=1}^{M-1} (e^{v_{\mu}} -1) P_{\mu} ]^{N},
\end{displaymath} (75)

erhalten wir
$\displaystyle <K_{\nu}>$ $\textstyle =$ $\displaystyle \frac{\partial M_{\vec{K}} (\vec{v}=0)}{\partial v_{\nu}}
= N P_{\nu},$ (76)
$\displaystyle <K_{\nu} K_{\mu}>$ $\textstyle =$ $\displaystyle \frac{\partial^{2} M_{\vec{K}}(\vec{v}=0)}
{\partial v_{\nu} \partial v_{\mu}} = N P_{\nu} \delta_{\nu\mu}+N (N-1) P_{\nu}
P_{\mu}$ (77)

und
\begin{displaymath}
<(K_{\nu}-<K_{\nu}>)(K_{\mu}-<K_{\mu}>)> = N P_{\nu} \delta_{\nu\mu}
- N P_{\nu} P_{\mu} .
\end{displaymath} (78)

Zur Vereinfachung der Schreibweise setzen wir
\begin{displaymath}
y_{\nu} = \frac{1}{\sqrt{N}} (K_{\nu} - N P_{\nu}).
\end{displaymath} (79)

Diese Größe besitzt die Mittelwerte und Kovarianzmatrix
$\displaystyle a_{\nu} \equiv <y_{\nu}>$ $\textstyle =$ $\displaystyle 0$ (80)
$\displaystyle (A^{-1})_{\nu\mu} \equiv <y_{\nu} y_{\mu} >$ $\textstyle =$ $\displaystyle P_{\nu} \delta_{\nu\mu}
- P_{\nu} P_{\mu}.$ (81)

Die Kovarianzmatrix $A^{-1}$ läßt sich leicht analytisch invertieren, und zwar
\begin{displaymath}
A_{\nu\mu} = \frac{1}{P_{M}} + \frac{\delta_{\nu\mu}}{P_{\mu}}.
\end{displaymath} (82)

Daraus erhalten wir die Verteilungsfunktion der Abweichungen $\vec{y}$ zu
\begin{displaymath}
p(y_{1},y_{2},...y_{M-1}) = C e^{-\frac{1}{2} \sum_{\mu=1}^{...
.../P_{\mu} - \frac{1}{2} (\sum_{\mu=1}^{M-1} y_{\mu})^{2}/P_{M}}
\end{displaymath} (83)

mit
\begin{displaymath}
C = \frac{1}{\sqrt{2^{M-1} \pi^{M-1} P_{1}P_{2}\cdot\cdot \cdot P_{M}}}.
\end{displaymath} (84)

In der obigen Ableitung haben wir benutzt, daß

\begin{displaymath}
y_{M} = - \sum_{\mu=1}^{M-1} y_{\mu}
\end{displaymath}

ist. Als neue Variable führen wir die Größe
\begin{displaymath}
\Theta \equiv z_{1} = \frac{1}{M-1} \sum_{\mu=1}^{M} y_{\mu}...
...^{2}/P_{\mu} +
(\sum_{\mu=1}^{M-1} y_{\mu})^{2}/P_{M} \right),
\end{displaymath} (85)

sowie die Hilfstransformationen
$\displaystyle z_{2}$ $\textstyle =$ $\displaystyle y_{2}$  
$\displaystyle z_{3}$ $\textstyle =$ $\displaystyle y_{3}$  
  $\textstyle .$    
  $\textstyle .$    
$\displaystyle z_{M-1}$ $\textstyle =$ $\displaystyle y_{M-1}$  

ein. Die Transformation der Verteilungsfunktion von der zufälligen Veränderlichen $\vec{y}$ auf die zufällige Veränderliche $\vec{z}$ und nachfolgende Integration (Reduktion) über die Hilfsgrößen $z_{2},z_{3},...,z_{M-1}$ führt nach einiger Rechnung auf die Dichtefunktion für die Veränderliche $\Theta \equiv z_{1}$:
\begin{displaymath}
p(\Theta) = \frac{1}{\tau^{n+1}} \Theta^{n} \frac{e^{-\Theta/\tau}}{n!},
\end{displaymath} (86)

mit

\begin{displaymath}
\tau= \frac{2}{M-1}, \; \; \; \; \; n=\frac{M-1}{2} -1 = \frac{1-\tau}{\tau}.
\end{displaymath}

Diese spezielle Gammafunktion nennt man auch die $\chi^{2}$- Verteilung. Sie hängt offensichtlich nur noch von der Anzahl $M$ der Histogramm- Intervalle ab, jedoch nicht von der ursprünglichen Dichtefunktion $p(\vec{a};x)$ bzw den Wahrscheinlichkeiten $P_{1},P_{2},...,P_{M}$. Die erzeugende Funktion der $\chi^{2}$- Verteilung ergibt sich zu
\begin{displaymath}
M_{\Theta}(v) = \frac{1}{(1-v\tau)^{n+1}}
= \frac{1}{(1-\frac{2v}{M-1})^{(M-1)/2-1}},
\end{displaymath} (87)

und die logarithmisch erzeugende Funktion zu
\begin{displaymath}
H_{\Theta}(v) = -(n+1) \; ln(1-vt) = -(\frac{M-1}{2} - 1) \; ln(1 -
\frac{2v}{M-1}).
\end{displaymath} (88)

Aus der Reihenentwicklung

\begin{displaymath}
H_{\Theta}(v) = (n+1) \sum_{\nu=1}^{\infty} \frac{\tau^{\nu}}{\nu} v^{\nu}
= v + \frac{1}{2} (\frac{2}{M-1}) v^{2} + .....
\end{displaymath}

erhalten wir die ersten beiden Momente der $\chi^{2}$- Verteilung:
$\displaystyle <\Theta>$ $\textstyle =$ $\displaystyle 1$ (89)
$\displaystyle <(\Theta - <\Theta>)^{2}>$ $\textstyle =$ $\displaystyle \frac{2}{M-1} .$ (90)

Dieses ist exakt das Ergebnis, das wir bereits im ersten Teil dieses Tutorials ohne Beweis benutzt hatten. Dort hatten wir eine weitere Größe eingeführt, nämlich

\begin{displaymath}
\Theta' = (\Theta -1) \sqrt{\tau} = (\Theta -1) \sqrt{\frac{2}{M-1}}.
\end{displaymath}

Die Dichtefunktion dieser Veränderlichen besitzt zwar die von $M$ unabhängigen ersten beiden Momente
$\displaystyle <\Theta'>$ $\textstyle =$ $\displaystyle 0$  
$\displaystyle <(\Theta' -<\Theta'>)^{2}>$ $\textstyle =$ $\displaystyle 1,$  

alle höheren Momente hängen aber weiterhin von $M$ ab. Da die Dichtefunktion der Veränderlichen $\Theta'$ wesentlich komplizierter als die der Größe $\Theta$ ist, benutzen wir im folgenden den $\chi^{2}$- Test mit der Veränderlichen $\Theta$. Die Größe $M-1$ nennt man die Anzahl der Freiheitsgrade (NDF = ''Number Degrees of Freedom'').

Der $\chi^{2}$- Test
Zusammenfassend geben wir noch einmal die Schritte an, die zum $\chi^{2}$- Test führen:
1. Aus einem Versuch erhalten wir $N$ unabhängige Messungen der zufälligen Veränderlichen $x$, und zwar $x_{1},x_{2},...,x_{N}$.
2. Wir sortieren die Messungen in ein Histogramm mit $M$ Intervallen zwischen den Stützstellen $x^{(\mu)}, \mu=0,1,2,...,M$. Die Häufigkeiten in den Intervallen seien $K_{1},K_{2},...,K_{M}$.
3. Aus der Hypothese, daß die Dichtefunktion $p(\vec{a};x)$ den Bedingungskomplex des Versuches vollständig beschreibt, erhalten wir die hypothetischen Wahrscheinlichkeiten

\begin{displaymath}
P_{\mu} = \int_{x^{(\mu-1)}}^{x^{(\mu)}} dx p(\vec{a};x).
\end{displaymath}

4. Wir bilden das sogenannte $\chi^{2}/NDF$,
\begin{displaymath}
\chi_{0}^{2}/NDF \equiv \Theta_{0} = \frac{1}{M-1} \sum_{\mu=1}^{M}
\frac{(K_{\mu}-N P_{\mu})^{2}}{N P_{\mu}},
\end{displaymath} (91)

und berechnen die Wahrscheinlichkeit dafür, daß die Veränderliche $\Theta$ einen Wert größer als den gemessenen Wert $\Theta_{0}$ annimmt:
\begin{displaymath}
P(\Theta > \Theta_{0}) = \int_{\Theta_{0}}^{\infty} d\Theta
\frac{1}{\tau^{n+1}} \Theta^{n} \frac{e^{-\Theta/\tau}}{n!},
\end{displaymath} (92)

mit

\begin{displaymath}
\tau = \frac{2}{M-1}, \; \; \; \; \; n=\frac{1-\tau}{\tau}.
\end{displaymath}

Diese Wahrscheinlichkeit, ausgedrückt in $\%$, nennt man gemeinhin auch den $CL$ (Confidence Level).
5. Das Ergebnis lautet: Mit einem Confidence Level von $CL$ $\%$ erfüllt die Dichtefunktion $p(\vec{a};x)$ den Bedingungskomplex des Versuches mit den Messungen $x_{1},x_{2},...,x_{N}$.

Der $\chi^{2}$- Test ist erweiterungsfähig auf den allgemeinen Fall, daß $N$ Messungen $y_{\nu}, \nu=1,2,...,M$, mit den Meßfehlern $d y_{\nu}$ bei den Koordinaten $x_{\nu}$ gegeben sind. Wir nehmen an, daß die Messungen durch eine Funktion der Form $y=f(x)$ beschrieben werden können und bilden das $\chi^{2}/NDF$,

\begin{displaymath}
\chi_{0}^{2}/NDF \equiv \Theta_{0} = \frac{1}{M-1} \sum_{\nu=1}^{M}
\frac{(y_{\nu}-f(x_{\nu}))^{2}}{(dy_{\nu})^{2}}.
\end{displaymath} (93)

Diese Veränderliche gehorcht allerdings nur dann der $\chi^{2}$- Verteilung, sofern die Meßfehler rein statistischer Natur sind.



Harm Fesefeldt
2006-05-09