diff options
Diffstat (limited to '')
-rw-r--r-- | buch/chapters/40-eigenwerte/eigenwerte.tex | 1 | ||||
-rw-r--r-- | buch/chapters/40-eigenwerte/normalformen.tex | 1 | ||||
-rw-r--r-- | buch/chapters/40-eigenwerte/spektralradius.tex | 5 | ||||
-rw-r--r-- | buch/chapters/60-gruppen/symmetrien.tex | 11 | ||||
-rw-r--r-- | buch/chapters/70-graphen/beschreibung.tex | 2 | ||||
-rw-r--r-- | buch/chapters/70-graphen/waerme.tex | 179 | ||||
-rw-r--r-- | buch/chapters/70-graphen/wavelets.tex | 55 | ||||
-rw-r--r-- | buch/chapters/80-wahrscheinlichkeit/chapter.tex | 8 | ||||
-rw-r--r-- | buch/chapters/80-wahrscheinlichkeit/google.tex | 155 | ||||
-rw-r--r-- | buch/chapters/80-wahrscheinlichkeit/markov.tex | 198 | ||||
-rw-r--r-- | buch/chapters/80-wahrscheinlichkeit/parrondo.tex | 74 | ||||
-rw-r--r-- | buch/chapters/80-wahrscheinlichkeit/positiv.tex | 185 | ||||
-rw-r--r-- | buch/chapters/90-crypto/chapter.tex | 10 | ||||
-rw-r--r-- | buch/chapters/references.bib | 32 | ||||
-rw-r--r-- | buch/papers/erdbeben/Teil_Fabio.tex | 238 | ||||
-rw-r--r-- | buch/papers/erdbeben/references.bib | 7 | ||||
-rw-r--r-- | buch/papers/erdbeben/teil0.tex | 138 | ||||
-rw-r--r-- | buch/papers/erdbeben/teil1.tex | 341 | ||||
-rw-r--r-- | buch/papers/munkres/teil1.tex | 2 |
19 files changed, 966 insertions, 676 deletions
diff --git a/buch/chapters/40-eigenwerte/eigenwerte.tex b/buch/chapters/40-eigenwerte/eigenwerte.tex index 1af91f8..f0d7b16 100644 --- a/buch/chapters/40-eigenwerte/eigenwerte.tex +++ b/buch/chapters/40-eigenwerte/eigenwerte.tex @@ -80,6 +80,7 @@ E_\lambda \{ v\;|\; Av=\lambda v\} \] der {\em Eigenraum} zum Eigenwert $\lambda$. +\index{Elambda(A)@$E_\lambda(A)$}% \index{Eigenraum}% \end{definition} diff --git a/buch/chapters/40-eigenwerte/normalformen.tex b/buch/chapters/40-eigenwerte/normalformen.tex index e59f1dc..96cb18b 100644 --- a/buch/chapters/40-eigenwerte/normalformen.tex +++ b/buch/chapters/40-eigenwerte/normalformen.tex @@ -103,6 +103,7 @@ ist mit $\lambda_i\in\Bbbk'$. Nach Satz~\ref{buch:eigenwerte:satz:zerlegung-in-eigenraeume} liefern die verallgemeinerten Eigenräume $V_i=\mathcal{E}_{\lambda_i}(A)$ eine +\index{Elambda@$\mathcal{E}_{\lambda}(A)$}% Zerlegung von $V$ in invariante Eigenräume \[ V=V_1\oplus V_2\oplus \dots\oplus V_l diff --git a/buch/chapters/40-eigenwerte/spektralradius.tex b/buch/chapters/40-eigenwerte/spektralradius.tex index 1cdaf35..c0d4de9 100644 --- a/buch/chapters/40-eigenwerte/spektralradius.tex +++ b/buch/chapters/40-eigenwerte/spektralradius.tex @@ -585,6 +585,7 @@ Dies führt uns auf die Grösse \limsup_{n\to\infty} \|M^n\|^\frac1n, \label{buch:eqn:gelfand-grenzwert} \end{equation} +\index{pi(M)@$\pi(M)$}% die darüber entscheidet, ob die Potenzreihe $f(A)$ konvergiert. @@ -631,9 +632,11 @@ Viel einfacher ist der Begriff des Spektralradius. \begin{definition} \label{buch:definition:spektralradius} -Der {\em Spektralradius} der Matrix $M$ ist der Betrag des betragsgrössten +Der {\em Spektralradius} $\varrho(M)$ der Matrix $M$ ist der Betrag des +betragsgrössten \index{Spektralradius}% Eigenwertes. +\index{rho(M)@$\varrho(M)$}% \end{definition} Wir wollen in diesem Abschnitt zeigen, dass der Gelfand-Radius mit diff --git a/buch/chapters/60-gruppen/symmetrien.tex b/buch/chapters/60-gruppen/symmetrien.tex index 252fdca..7222c2c 100644 --- a/buch/chapters/60-gruppen/symmetrien.tex +++ b/buch/chapters/60-gruppen/symmetrien.tex @@ -4,6 +4,7 @@ % % (c) 2020 Prof Dr Andreas Müller, Hochschule Rapperswil % +% !TeX spellcheck = de_CH \section{Symmetrien \label{buch:section:symmetrien}} \rhead{Symmetrien} @@ -21,7 +22,7 @@ Begriffs verständlich macht. \centering \includegraphics[width=\textwidth]{chapters/60-gruppen/images/castle.jpeg} \caption{Das Castle Howard in Yorkshire war in dieser ausgeprägt symmetrischen -Form geplant, wurde dann aber in modifizeirter Form gebaut. +Form geplant, wurde dann aber in modifizierter Form gebaut. Messungen zwischen Punkten in der rechten Hälfte des Bildes ergeben die gleichen Werte wie Messungen entsprechenden Strecken in der linken Hälfte, was den Begriff Symmetrie rechtfertigt. @@ -222,7 +223,7 @@ x(t) = \frac{1}{\omega} \sin\omega t, \qquad p(t) = \cos \omega t. \] -In Matrixform kann man die allgemeine Lösung zur Anfangsbedingun $x(0)=x_0$ +In Matrixform kann man die allgemeine Lösung zur Anfangsbedingung $x(0)=x_0$ und $p(0)=p_0$ \begin{equation} \begin{pmatrix} @@ -274,7 +275,7 @@ $\operatorname{GL}_n(\mathbb{R})$, da gilt. Die Lösungen der Differentialgleichung~\eqref{chapter:gruppen:eqn:phasenraumdgl} -sind in Abbildung~\ref{chapter:gruppen:fig:phasenraum} +sind in Abbildung~\ref{chapter:gruppen:fig:phasenraum} dargestellt. dargestellt. Die Matrizen $\Phi_t$ beschreiben eine kontinuierliche Symmetrie des Differentialgleichungssystems, welches den harmonischen Oszillator @@ -408,7 +409,7 @@ Auch die Ableitung funktioniert dort nicht mehr. Bewegt man sich mit konstanter Geschwindigkeit über den Nordpol, springt die Ableitung der geographischen Breite von einem positiven Wert auf einen negativen Wert, sie kann also nicht differenzierbar sein. -Diese Einschränkungen sind in der Praxis nur ein geringes Problem dar, +Diese Einschränkungen stellen in der Praxis nur ein geringes Problem dar, da die meisten Reisen nicht über die Pole erfolgen. Der Polarforscher, der in unmittelbarer Umgebung des Poles arbeitet, @@ -661,7 +662,7 @@ $\gamma(t)$ ist $\dot{y}(t)$ in den Karten $\varphi_1$ und $\varphi_2$ und $\dot{x}(t)$ in den Karten $\varphi_3$ und $\varphi_4$. Die spezielle Kurve $\gamma(t) = (\cos t,\sin t)$ hat in einem Punkt -$t\in (0,\frac{\pi}2)$. +$t\in (0,\frac{\pi}2)$ in der Karte $\varphi_1$ den Tangentialvektor $\dot{y}(t)=\cos t$, in der Karte $\varphi_3$ aber den Tangentialvektor $\dot{x}=-\sin t$. Die Ableitung des Kartenwechsels in diesem Punkt ist die $1\times 1$-Matrix diff --git a/buch/chapters/70-graphen/beschreibung.tex b/buch/chapters/70-graphen/beschreibung.tex index 918594d..af934e4 100644 --- a/buch/chapters/70-graphen/beschreibung.tex +++ b/buch/chapters/70-graphen/beschreibung.tex @@ -193,7 +193,7 @@ a_{i\!j} 1&\qquad (j,i) \in E\\ 0&\qquad \text{sonst.} \end{cases} -\label{buch:graphen:eqn:adjazenzmatrix} +\label{buch:graphen:eqn:adjazenzmatrixgerichtet} \end{equation} Die Matrix $A(G)$ hat also genau dann einen nicht verschwindenden Matrixeintrag in Zeile $i$ und Spalte $j$, wenn es eine Verbindung diff --git a/buch/chapters/70-graphen/waerme.tex b/buch/chapters/70-graphen/waerme.tex index e7fc023..bfeff74 100644 --- a/buch/chapters/70-graphen/waerme.tex +++ b/buch/chapters/70-graphen/waerme.tex @@ -5,10 +5,11 @@ % \section{Wärmeleitung auf einem Graphen \label{buch:section:waermeleitung-auf-einem-graphen}} -Die Vektoren, auf denen die Laplace-Matrix operiert, können betrachtet -werden als Funktionen, die jedem Knoten einen Wert zuordnen. +Die Vektoren, auf denen die Laplace-Matrix operiert, können +als Funktionen betrachtet werden, die jedem Knoten einen Wert zuordnen. Eine mögliche physikalische Interpretation davon ist die Temperaturverteilung auf dem Graphen. +\index{Temperaturverteilung}% Die Kanten zwischen den Knoten erlauben der Wärmeenergie, von einem Knoten zu einem anderen zu fliessen. Je grösser die Temperaturdifferenz zwischen zwei Knoten ist, desto @@ -29,7 +30,7 @@ d_iT_i \biggr) \] Der Term auf der rechten Seite ist genau die Wirkung der -Laplace-Matrix auf dem Vektor $T$ der Temperaturen: +Laplace-Matrix $L=L(G)$ auf dem Vektor $T$ der Temperaturen: \begin{equation} \frac{dT}{dt} = @@ -38,6 +39,7 @@ Laplace-Matrix auf dem Vektor $T$ der Temperaturen: \end{equation} Der Wärmefluss, der durch die Wärmeleitungsgleichung~\eqref{buch:graphen:eqn:waermeleitung} beschrieben +\index{Wärmeleitungsgleichung}% wird, codiert ebenfalls wesentliche Informationen über den Graphen. Je mehr Kanten es zwischen verschiedenen Teilen eines Graphen gibt, desto schneller findet der Wärmeaustausch zwischen diesen Teilen @@ -50,6 +52,7 @@ Die Lösungen der Wärmeleitungsgleichung liefern also Informationen Die Wärmeleitungsgleichung~\eqref{buch:graphen:eqn:waermeleitung} ist eine lineare Differentialgleichung mit konstanten Koeffizienten, die mit der Matrixexponentialfunktion gelöst werden. +\index{Matrixexponentialfunktion}% Die Lösung ist \[ f(t) = e^{-\kappa Lt}f(0). @@ -58,93 +61,131 @@ f(t) = e^{-\kappa Lt}f(0). Die Berechnung der Lösung mit der Matrixexponentialreihe ist ziemlich ineffizient, da grosse Matrizenprodukte berechnet werden müssen. Da die Matrix $L$ symmetrisch ist, gibt es eine Basis aus -orthonormierten Eigenvektoren und die Eigenwerte sind reell. -Wir bezeichnen die Eigenvektoren mit $f_1,\dots,f_n$ und die +orthonormierten Eigenvektoren und die zugehörigen Eigenwerte sind reell. +Wir bezeichnen die Eigenvektoren mit $\chi_1,\dots,\chi_n$ und die zugehörigen Eigenwerte mit $\lambda_i$. -Die Funktion $f_i(t)= e^{-\kappa\lambda_it}f_i$ ist dann eine Lösung +Die Funktion $\chi_i(t)= e^{-\kappa\lambda_it}\chi_i$ ist dann eine Lösung der Wärmeleitungsgleichung, denn die beiden Seiten -\begin{align*} -\frac{d}{dt}f_i(t) +\begin{equation} +\begin{aligned} +\text{linke Seite:}&& +\frac{d}{dt}\chi_i(t) &= --\kappa\lambda_ie^{-\kappa\lambda_it}f_i +-\kappa\lambda_ie^{-\kappa\lambda_it}\chi_i = --\kappa\lambda_i f_i(t) +-\kappa\lambda_i \chi_i(t) \\ --\kappa Lf_i(t) +\text{rechte Seite:}&& +-\kappa L\chi_i(t) &= --\kappa e^{-\kappa\lambda_it} Lf_i +-\kappa e^{-\kappa\lambda_it} L\chi_i = --\kappa e^{-\kappa\lambda_it} \lambda_i f_i +-\kappa e^{-\kappa\lambda_it} \lambda_i \chi_i = --\kappa \lambda_i f_i(t) -\end{align*} +-\kappa \lambda_i \chi_i(t) +\end{aligned} +\end{equation} von \eqref{buch:graphen:eqn:waermeleitung} stimmen überein. Eine Lösung der Wärmeleitungsgleichung zu einer beliebigen Anfangstemperaturverteilung $f$ kann durch Linearkombination aus -den Lösungen $f_i(t)$ zusammengesetzt werden. -Dazu ist nötig, $f$ aus den Vektoren $f_i$ linear zu kombinieren. -Da aber die $f_i$ orthonormiert sind, ist dies besonders einfach, +den Lösungen $\chi_i(t)$ zusammengesetzt werden. +Dazu ist nötig, $f$ aus den Vektoren $\chi_i$ linear zu kombinieren. +Da aber die $\chi_i$ orthonormiert sind, ist dies besonders einfach, die Koeffizienten sind die Skalarprodukte mit den Eigenvektoren: \[ -f=\sum_{i=1}^n \langle f_i,f\rangle f_i. +f=\sum_{i=1}^n \langle \chi_i,f\rangle \chi_i. \] -Daraus kann man die allgmeine Lösungsformel +Daraus kann man die allgemeine Lösungsformel \begin{equation} f(t) = -\sum_{i=1}^n \langle f_i,f\rangle f_i(t) +\sum_{i=1}^n \langle \chi_i,f\rangle \chi_i(t) = -\sum_{i=1}^n \langle f_i,f\rangle e^{-\kappa\lambda_i t}f_i +\sum_{i=1}^n \langle \chi_i,f\rangle e^{-\kappa\lambda_i t}\chi_i \label{buch:graphen:eqn:eigloesung} \end{equation} ableiten. -\subsection{Beispiel: Ein zyklischer Graph} +\subsection{Beispiel: Ein zyklischer Graph +\label{buch:graphen:subsection:zyklischer-graph}} \begin{figure} \centering \includegraphics{chapters/70-graphen/images/kreis.pdf} -\caption{Beispiel Graph zur Illustration der verschiedenen Basen auf einem +\caption{Beispielgraph zur Illustration der verschiedenen Basen auf einem Graphen. \label{buch:graphen:fig:kreis}} \end{figure} Wir illustrieren die im folgenden entwickelte Theorie an dem Beispielgraphen von Abbildung~\ref{buch:graphen:fig:kreis}. -Besonders interessant sind die folgenden Funktionen: +Für jedes $k=0,\dots,n-1$ ist der Vektor mit den Komponenten \[ -\left. -\begin{aligned} -s_m(k) +\chi_k(l) = e^{2\pi ikl/n}, \quad l=1,\dots,n +\] +ein Eigenvektor der Laplace-Matrix zum Eigenwert +$\lambda_k=4\sin^2\frac{\pi k}{n}$. +Tatsächlich ist +\begin{align*} +(L\chi_k)(l) &= -\sin\frac{2\pi mk}{n} +-\chi_k(l-1) ++ +2\chi_k(l) +- +\chi_k(l+1) \\ -c_m(k) &= -\cos\frac{2\pi mk}{n} -\end{aligned} -\; -\right\} -\quad -\Rightarrow -\quad -e_m(k) +-e^{2\pi ik(l-1)/n} ++ +2e^{2\pi ikl/n} +- +e^{2\pi ik(l+1)/n} +\\ +&= +(-e^{-2\pi ik/n}+2-e^{2\pi ik/n})e^{2\pi ikl/n} +\\ +&= +-(e^{2\pi ik/2n}-e^{-2\pi ik/2n})^2 \chi_k(l) +\\ +&= +- +\biggl( +\frac{e^{2\pi ik/2n}-e^{-2\pi ik/2n}}{2i} +\biggr)^2 +(2i)^2 \chi_k(l) +\\ +&= +4\sin^2\frac{\pi k}n \chi_k(l) +\end{align*} + +Natürlich sind auch Real- und Imaginärteil Eigenvektoren: +\[ +\begin{aligned} +s_k(l) +&= +\sin\frac{2\pi kl}{n} = -e^{2\pi imk/n} +\Im \chi_k(l) +\\ +c_k(l) +&= +\cos\frac{2\pi kl}{n} = -c_m(k) + is_m(k). +\Re\chi_k(l) +\end{aligned} \] Das Skalarprodukt dieser Funktionen ist \[ -\langle e_m, e_{m'}\rangle +\langle \chi_m, \chi_{m'}\rangle = \frac1n -\sum_{k=1}^n -\overline{e^{2\pi i km/n}} -e^{2\pi ikm'/n} +\sum_{l=1}^n +\overline{e^{2\pi i ml/n}} +e^{2\pi im'l/n} = \frac1n -\sum_{k=1}^n -e^{\frac{2\pi i}{n}(m'-m)k} +\sum_{l=1}^n +e^{\frac{2\pi i}{n}(m'-m)l} = \delta_{mm'} \] @@ -157,21 +198,9 @@ c_0, c_1,s_1,c_2,s_2,\dots c_{\frac{n}2-1},c_{\frac{n}2-1},c_{\frac{n}2} \] eine orthonormierte Basis. - -Die Laplace-Matrix kann mit der folgenden Definition zu einer linearen -Abbildung auf Funktionen auf dem Graphen gemacht werden. -Sei $f\colon V\to \mathbb{R}$ und $L$ die Laplace-Matrix mit -Matrixelementen $l_{vv'}$ wobei $v,v'\in V$ ist. -Dann definieren wir die Funktion $Lf$ durch -\[ -(Lf)(v) -= -\sum_{v'\in V} l_{vv'}f(v'). -\] - \subsection{Standardbasis und Eigenbasis \label{buch:subsection:standardbasis-und-eigenbasis}} -Die einfachste Basis, aus der siche Funktionen auf dem Graphen linear +Die einfachste Basis, aus der sich Funktionen auf dem Graphen linear kombinieren lassen, ist die Standardbasis. Sie hat für jeden Knoten $v$ des Graphen eine Basisfunktion mit den Werten \[ @@ -180,5 +209,37 @@ e_v\colon V\to\mathbb R:v'\mapsto \begin{cases} 0\qquad&\text{sonst.} \end{cases} \] +Sie zeichnet sich dadurch aus, dass sie perfekt lokalisiert ist. +Im Gegensatz dazu zeigt das Beispiel von +Abschnitt~\ref{buch:graphen:subsection:zyklischer-graph}, dass +die Eigenfunktionen von $L(G)$ typischerweise delokalisiert sind. +Im Beispiel hat $\chi_k(l)$ überall auf dem Graphen den gleichen +Betrag. +Die ``Frequenz'' einer Eigenfunktion dagegen ist exakt bestimmt. + +\subsection{Fourier-Theorie auf einem Graphen} +Die Eigenfunktionen der Laplace-Matrix auf einem Graphen erlauben +also, das Wärmeleitungsproblem auf dem Graphen auf ganz ähnliche +Art zu lösen, wie die Fourier-Theorie das Wärmeleitungsproblem auf +$\mathbb{R}$ oder auf einem Intervall löst. +Es ist daher angemessen, die Entwicklung einer Funktion +$f\colon G\to\mathbb{C}$ nach den Eigenvektoren $\chi_k$ +als Fourier-Transformation zu bezeichnen und die Koeffizienten +\( +c_k = \langle \chi_k, f\rangle +\) +als die Fourier-Koeffizienten. +Grundlegende Eigenschaften der Fourier-Transformation stehen damit +auch für die Analyse von Funktionen auf einem Graphen zur Verfügung. +Es fehlen allerdings Eigenschaften, die mit zusätzlicher Struktur +auf dem Definitionsbereich zusammenhängen. +Die Faltung zum Beispiel setzt eine Rechenoperation auf dem +Definitionsbereich voraus, welche natürlich in einem Graphen nicht erwartet +werden kann. +Im Beispiel von Abschnitt~\ref{buch:graphen:subsection:zyklischer-graph} +lässt sich eine solche Struktur finden, die Knoten des Graphen können +als die Elemente einer zyklischen Gruppe betrachtet werden. +Daraus lassen sich die bekannten Faltungsformeln der diskreten +Fourier-Transformation ableiten. diff --git a/buch/chapters/70-graphen/wavelets.tex b/buch/chapters/70-graphen/wavelets.tex index 2b9f29b..b982bce 100644 --- a/buch/chapters/70-graphen/wavelets.tex +++ b/buch/chapters/70-graphen/wavelets.tex @@ -60,7 +60,7 @@ Lösungenfunktionen durch Linearkombination erzeugt werden können. In der Standarbasis (links) ist es am einfachsten, die Funktionswerte abzulesen, in der Eigenbasis (Mitte) kann die zeitliche Entwicklung besonders leicht berechnet werden. -Dazuwischen liegen die Fundamentallösungen (rechts), die eine einigermassen +Dazwischen liegen die Fundamentallösungen (rechts), die eine einigermassen übersichtliche Zeitentwicklung haben, die Berechnung der Temperatur an einer Stelle $x$ zur Zeit $t$ ist aber erst durch das Integral \eqref{buch:graphen:eqn:fundamentalueberlagerung} gegeben. @@ -73,20 +73,21 @@ Standardbasisvektor mit Hilfe der Lösungsformel~\eqref{buch:graphen:eqn:eigloesung} gefunden werden. Aus physikalischen Gründen ist aber offensichtlich, dass die -Wärmeenergie Fundamentallösungen $F_i(t)$ für kurze Zeiten $t$ -in der Nähe des Knoten $i$ konzentriert ist. -Dies ist aber aus der expliziten Formel +Wärmeenergie der Fundamentallösungen $F_i(t)$ für kurze Zeiten $t$ +in der Nähe des Knotens $i$ konzentriert ist. +Dies ist aber aus der Fourier-Entwicklung \begin{equation} F_i(t) = -\sum_{j=1}^n \langle f_j,e_i\rangle e^{-\kappa \lambda_i t} f_j +\sum_{j=1}^n \langle \chi_j,e_i\rangle e^{-\kappa \lambda_i t} \chi_j = \sum_{j=1}^n \overline{f}_{ji} e^{-\kappa \lambda_i t}, \label{buch:graphen:eqn:fundamentalgraph} \end{equation} nicht unmittelbar erkennbar. -Man kann aber aus~\eqref{buch:graphen:eqn:fundamentalgraph} ablesen, +Man kann aber aus~\eqref{buch:graphen:eqn:fundamentalgraph} +wenigstens ablesen, dass für zunehmende Zeit die hohen Frequenzen sehr schnell gedämpft werden. Die hohen Frequenzen erzeugen also den scharfen Peak für Zeiten nahe @@ -115,7 +116,7 @@ Die Darstellung im Frequenzraum und in der Zeit sind also extreme Darstellungen, entweder Frequenzlokalisierung oder zeitliche Lokalisierung ermöglichen, sich aber gegenseitig ausschliessen. -\subsubsection{Dilatation} +\subsubsection{Dilatation im Frequenzraum, spektrale Dilatation} Eine Wavelet-Basis für die $L^2$-Funktionen auf $\mathbb{R}$ erlaubt eine Funktion auf $\mathbb{R}$ auf eine Art zu analysieren, die eine ungenaue zeitliche Lokalisierung bei entsprechend ungenauer @@ -140,7 +141,7 @@ Graphen gibt es keine Rechtfertigung für diese spezielle Wahl von Streckungsfaktoren mehr. Es stellt sich daher die Frage, ob man für eine beliebige Menge \( -T= \{ t_1,t_2,\dots\} \} +T= \{ t_1,t_2,\dots\} \) von Streckungsfaktoren eine Familie von Funktionen $\chi_j$ zu finden derart, dass man sich die $\chi_j$ in einem gewissen Sinn als aus @@ -164,14 +165,14 @@ Menge von reellen Zahlen ohne innere algebraische Struktur ist. \centering \includegraphics{chapters/70-graphen/images/gh.pdf} \caption{Lokalisierungsfunktion $g(\lambda)$ für die Dilatation (links). -Die Dilatierten Funktionen $g_i=\tilde{D}_{1/a_i}g$ lokalisieren +Die dilatierten Funktionen $g_i=\tilde{D}_{1/a_i}g$ lokalisieren die Frequenzen jeweils um die Frequenzen $a_i$ im Frequenzraum. Der Konstante Vektor ist vollständig delokalisiert, die Funktion $h$ in der rechten Abbildung entfernt die hohen Frequenzen und liefert Funktionen, -die in der Umgebung eines Knotens wie die Konstante Funktion aussehen. +die in der Umgebung eines Knotens wie die konstante Funktion aussehen. \label{buch:graphs:fig:lokalisierung}} \end{figure} -Das Mutter-Wavelet einer Wavelet-Analyse zeichnet definiert, in welchem Mass +Das Mutter-Wavelet einer Wavelet-Analyse definiert, in welchem Mass sich Funktionen im Orts- und im Frequenzraum lokalisieren lassen. Die Standardbasis der Funktionen auf einem Graphen repräsentieren die perfekte örtliche Lokalisierung, Eigenbasis der Laplace-Matrix $L$ repräsentiert @@ -181,8 +182,8 @@ $\lambda\to\infty$ rasch abfällt mit einem Maximum irgendwo dazwischen (Abbildung~\ref{buch:graphs:fig:lokalisierung}). Sie kann als eine Lokalisierungsfunktion im Frequenzraum betrachtet werden. -Die Matrix $g(L)$ bildet entfernt aus einer Funktion die ganz hohen und -die ganz tiefen Frequenz, lokalisiert also die Funktionen im Frequenzraum. +Die Matrix $g(L)$ entfernt die ganz hohen und die ganz tiefen Frequenz +aus einer Funktion, lokalisiert also die Funktionen im Frequenzraum. Die Standardbasisvektoren werden dabei zu Funktionen, die nicht mehr nur auf einem Knoten von $0$ verschieden sind, aber immer noch einigermassen auf dem Graphen lokalisiert sind. @@ -191,7 +192,7 @@ $\lambda_0 < \lambda_1\le \dots\le \lambda_n$ der Laplace-Matrix von Interesse. Die Matrix $g(L)$ kann mit Hilfe der Spektraltheorie berechnet werden, -was im vorliegenden Fall naheliegend ist, weil ja die Eigenvektoren von +was im vorliegenden Fall naheliegend ist, weil ja die Eigenvektoren der Laplace-Matrix bereits bekannt sind. Die Matrix $\chi^t$ bildet die Standardbasisvektoren in die Eigenbasis-Vektoren ab, also in eine Zerlegung im Frequenzraum ab, @@ -211,7 +212,7 @@ g(\lambda_0)&0&\dots&0\\ \label{buch:graphen:eqn:mutterwavelet} \end{equation} -\subsubsection{Dilatation} +\subsubsection{Spektrale Dilatation der Mutterwavelets} Die Dilatation um $a$ im Ortsraum wird zu einer Dilatation um $1/a$ im Frequenzraum. Statt also nach einer echten Dilatation der Spaltenvektoren in $g(L)$ @@ -266,12 +267,20 @@ h(L) + \sum_{i}g_i(L)=I gelten würde. Nach der Spektraltheorie gilt das nur, wenn für alle Eigenwerte $\lambda_k$, $k=1,\dots,n$ -\[ +\begin{equation} h(\lambda_k) + \sum_ig(a_i\lambda_k)=1 -\] +\label{buch:graphen:eqn:summegh} +\end{equation} gilt. -Für beliebige Funktionen $g$ und $h$ kann man nicht davon ausgehen, -aber man kann erwarten. + +Allerdings kann man im Allgemeinen nicht erwarten, +dass \ref{buch:graphen:eqn:summegh} für +beliebige Funktionen $g$ und $h$ gilt. +Da es aber nur auf die Werte auf den Eigenwerten ankommt, +muss nur sichergestellt sein, dass +die linke Seite von \eqref{buch:graphen:eqn:summegh} +nicht verschwindet. +Dies garantiert, dass die Wavelet-Entwicklung umkehrbar ist. Man muss daher zusätzlich verlangen, dass \[ h(\lambda_k) + \sum_{i} g(a_i\lambda_k) > 0 @@ -301,7 +310,7 @@ B\|v\|^2 Die Zahlen $A$ und $B$ heissen die {\em Frame-Konstanten} des Frames. \end{definition} -Die oben gefundenen Vektoren, die Spalten Vektoren von $h(L)$ und $g_i(L)$ +Die oben gefundenen Vektoren, die Spaltenvektoren von $h(L)$ und $g_i(L)$, bilden daher ein Frame. Die Frame-Konstanten kann man unmittelbar ausrechnen. Der mittlere Term von \eqref{buch:graphen:eqn:frame} ist @@ -318,12 +327,14 @@ h(\lambda)^2 + \sum_i g_i(\lambda)^2 \] abgeschätzt werden kann. Die Frame-Konstanten sind daher -\begin{align*} +\[ +\begin{aligned} A&=\min_{k} f(\lambda_k) & &\text{und}& B&=\max_{k} f(\lambda_k). -\end{align*} +\end{aligned} +\] Die Konstruktion hat also ein Frame für die Funktionen auf dem Graphen etabliert, die viele Eigenschaften einer Multiskalenanalyse in diese wesentlich weniger symmetrische Situation rettet. diff --git a/buch/chapters/80-wahrscheinlichkeit/chapter.tex b/buch/chapters/80-wahrscheinlichkeit/chapter.tex index 85b6d8c..270c44a 100644 --- a/buch/chapters/80-wahrscheinlichkeit/chapter.tex +++ b/buch/chapters/80-wahrscheinlichkeit/chapter.tex @@ -12,8 +12,8 @@ jedem Vektor einen neuen Vektor zuordnet. Es ist daher nicht abwegig zu erwarten, dass sich die Zeitentwicklung eines vom Zufall beeinflussten Systems, welches sich in mehreren verschiedenen Zuständen befinden kann, ebenfalls mit Hilfe -von Matrizen beschreiben lässt. -Eine solche Beschreiben ermöglicht leicht Verteilungen, +von Matrizen modellieren lässt. +Eine solche Beschreibung ermöglicht Verteilungen, Erwartungswerte und stationäre Zustände zu ermitteln. Im Abschnitt~\ref{buch:section:google-matrix} wird an Hand der Google @@ -22,6 +22,10 @@ auf eine Matrix führt. Abschnitt~\ref{buch:section:diskrete-markov-ketten} stellt dann die abstrakte mathematische Theorie der Markov-Ketten dar und behandelt einige wichtige Eigenschaften von Wahrscheinlichkeitsmatrizen. +Sie folgen aus den Eigenschaften positiver Matrizen und Vektoren, +die die Theorie von Perron und Frobenius in +Abschnitt~\ref{buch:section:positive-vektoren-und-matrizen} +beweist. Es stellt sich heraus, dass thermodynamische Quantensysteme sehr gut mit solchen Matrizen beschrieben werden können, zum Beispiel kann man einfache Formen von Laser auf diese Art behandeln. diff --git a/buch/chapters/80-wahrscheinlichkeit/google.tex b/buch/chapters/80-wahrscheinlichkeit/google.tex index ca78b3d..c9d0d8c 100644 --- a/buch/chapters/80-wahrscheinlichkeit/google.tex +++ b/buch/chapters/80-wahrscheinlichkeit/google.tex @@ -8,6 +8,7 @@ \rhead{Google-Matrix} Das Internet besteht aus einer grossen Zahl von Websites, etwa 400~Millionen aktiven Websites, jede besteht aus vielen einzelnen Seiten. +\index{Internet}% Es ist daher angemessen von $N\approx 10^9$ verschiedenen Seiten auszugehen. Eine natürliche Sprache umfasst dagegen nur einige 100000 bis Millionen von Wörtern. @@ -17,21 +18,23 @@ Die Zahl der Wörter, die auf den $N$ Seiten vorkommen können, ist also viel kleiner als die Zahl der zur Verfügung stehenden Wörter. Ein einzelnes Wort wird daher notwendigerweise auf einer grossen Zahl von Seiten vorkommen. -Eine Suche nach einem bestimmten Wort wird also in der überwiegenden Zahl +Eine Suche nach einem bestimmten Wort wird in der überwiegenden Zahl der Fälle derart viele Treffer zurückgeben, dass das Suchresultat nur dann nützlich sein kann, wenn eine zusätzliche Informationsquelle -ermöglicht, die Treffer in eine sinnvolle Ordnung zu bringem. +ermöglicht, die Treffer in eine sinnvolle Ordnung zu bringen. Genau dieses Problem stellte sich den vielen traditionellen Suchmaschienen in der ersten grossen Boomphase des Internets. -Traditionelle Informatione-Retrieval-Systeme operieren auf einem relativ +Traditionelle Information-Retrieval-Systeme operieren auf einem relativ +\index{Information-Retrieval}% kleinen Dokumentbestand und gehen davon aus, dass bereits wenige, spezifische Wörter nur in einem kleinen Teil des Dokumentbestandes vorkommen und damit eine übersichtliche Treffermenge ergeben. -Die Einengung der Treffermenge dank der Suche nach spezifischer Menge +Die Einengung der Treffermenge dank der Suche nach einzelnen Wörtern bedeutet aber auch, dass nach Synonymen oder alternative Formen eines Wortes separat gesucht werden muss, was die Übersichtlichkeit wieder zerstört. +\index{Treffermenge}% % % Ein Modell für Webseitenbesucher @@ -45,34 +48,35 @@ zerstört. \label{buch:figure:modellinternet}} \end{figure} -Das kombinierte Vorkommen von Wörtern oder Begriffen alleine kann also -nicht ausreichen, um die Seiten zum Beispiel einem Fachgebiet zuzuordnen. +Selbst das kombinierte Vorkommen von Wörtern oder Begriffen alleine reicht +nicht aus, um die Seiten zum Beispiel einem Fachgebiet zuzuordnen. Dazu muss eine externe Informationsquelle angezapft werden. Bei traditionellen Dokumenten liefert der Kontext, in dem ein Dokument erfasst wurde, solche ergänzenden Informationen. Eine Publikation in einem Fachjournal ordnet einen Text einem Fachgebiet zu. Im World-Wide-Web liefert die Link-Struktur diesen Kontext. -Dokumente zu ähnlichen Themen werden bevorzugt untereinander verlinkt -sein. +\index{Link}% +Dokumente zu ähnlichen oder verwandten Themen werden bevorzugt +untereinander verlinkt sein. Gesucht ist jetzt also ein Modell, welches objektiv die Linkstruktur -bewertet und daraus eine Rangordnung der passenden Wörter ableitet. +bewertet und daraus eine Rangordnung der Suchresultate ableitet. Die Linkstruktur kann natürlich als gerichteter Graph betrachtet und -mit Hilfe der Matrix~\eqref{buch:graphen:eqn:linkmatrix} -beschrieben werden. +mit Hilfe der Adjazenzmatrix~\eqref{buch:graphen:eqn:adjazenzmatrixgerichtet} +\index{Adjazenzmatrix}% +eines gerichteten Graphen beschrieben werden. Dies trägt jedoch der Anzahl der Wahlmöglichkeiten nicht Rechnung. -Eine Website mit nur einem Link auf die Seite $j$ hat mehr Gewicht -als eine Seite mit vielen Links, unter denen der Link auf die Seite $j$ -einer von vielen ist. -Im Beispiel-Inter der Abbildung~\ref{buch:figure:modellinternet} -signalisiert die Seite $t$ mit nur einem Link auf die Seite $8$ +Eine Website mit nur einem Link auf die Seite $j$ gibt der Seite $j$ +mehr Gewicht als eine Seite mit vielen Links, unter denen der Link +auf die Seite $j$ einer von Vielen ist. +Im Beispiel-Internet der Abbildung~\ref{buch:figure:modellinternet} +signalisiert die Seite $6$ mit nur einem Link auf die Seite $8$ viel deutlicher, dass $8$ eine wichtige Seite ist, also die die Seite $5$ tut, die auch noch zwei andere Links enthält. Wir können diesen Unterschied berücksichtigen, indem wir zu einem Wahrscheinlichkeitsmodell übergehen, was wir im folgenden Abschnitt tun werden. - % % Wahrscheinlichkeitsinterpretation % @@ -104,7 +108,9 @@ lange Zeit ungefähr nach den gleichen Dingen suchen und sich daher auf die gleiche Weise auf den verschiedenen Seiten verteilen und dass insbesondere die Verteilung stationär ist, dass also $P(S_i) = P(S'_i)$ gilt. +\index{Suchmaschine}% Suchmaschinen wie Google gehen davon aus, dass alle Besucher ungefähr +\index{Google}% die gleichen Suchprioritäten haben, so dass es sich lohnt, die Suchresultate nach der Wahrscheinlichkeit $P(S_i)$ zu ordnen und dem Suchenden die wahrscheinlichsten Dokumente als erste zu zeigen. @@ -113,19 +119,19 @@ wahrscheinlichsten Dokumente als erste zu zeigen. Um einen Zusammenhang zwischen $P(S_i)$ und $P(S'_j)$ herzustellen, muss die Navigation entlang der Links modelliert werden. Die naheliegende Wahrscheinlichkeitsinterpretation ist die bedingte -Wahrscheinlichkeit $P(S'_j|S_i)$ dass der Besucher auf der Seite $j$ +Wahrscheinlichkeit $P(S'_j\mid S_i)$ dass der Besucher auf der Seite $j$ landet, nachdem er auf der Seite $i$ die Linknavigation verwendet hat. Wenn es keinen Link zwischen den Seiten $i$ und $j$ gibt, dann ist diese -Navigation natürlich nicht möglich und es folgt $P(S'_j|S_i)=0$. -Falls es einen Link gibt, ist $P(S'_j|S_i)\ge 0$. +Navigation natürlich nicht möglich und es folgt $P(S'_j\mid S_i)=0$. +Falls es einen Link gibt, ist $P(S'_j\mid S_i)\ge 0$. A priori wissen wir nicht, wie wahrscheinlich es ist, dass der Besucher dem Link auf die Seite $j$ folgt, normalerweise werden nicht alle Links mit gleicher Wahrscheinlichkeit verwendet. -Wir nehmen daher zusätzlich an, dass alle Links gleich wahrscheinlich +Wir nehmen daher vereinfachend an, dass alle Links gleich wahrscheinlich sind. -Die Seite $i$ enthält $n_i$ Links, also ist die Wahrscheinlichkeit, -auf einer von $i$ aus verlinkten Seite $j$ zu landen $P(S'_j|S_i) = 1/n_i$. +Enthält die Seite $i$ genau $n_i$ Links, dann ist die Wahrscheinlichkeit, +auf einer von $i$ aus verlinkten Seite $j$ zu landen, $P(S'_j\mid S_i) = 1/n_i$. \subsubsection{Totale Wahrscheinlichkeit} Der Satz von der totalen Wahrscheinlichkeit ermöglicht, einen Zusammenhang @@ -136,13 +142,16 @@ Es gilt \begin{equation} P(S'_j) = -P(S'j|S_1) P(S_1) +P(S'j\mid S_1) P(S_1) + -P(S'j|S_2) P(S_2) +P(S'j\mid S_2) P(S_2) + \dots + -P(S'j|S_N) P(S_N). +P(S'j\mid S_N) P(S_N) += +\sum_{i=1}^N P(S_j'\mid S_i)P(S_i) +. \label{buch:google:eqn:totalewahrscheinlichkeit} \end{equation} Dies kann in Matrix- und Vektorform übersichtlicher geschrieben werden. @@ -170,18 +179,21 @@ P(S'_N) \end{pmatrix} \] zusammen. -Die bedingten Wahrscheinlichkeiten $h_{ji}=P(S'_j|S_i)$ sind mit zwei Indizes -beschrieben, sie bilden daher in natürlicher Weise eine Matrix -\[ +Die bedingten Wahrscheinlichkeiten $h_{ji}=P(S'_j\mid S_i)$ sind mit zwei Indizes +beschrieben, sie bilden daher in natürlicher Weise die sogenannte +{\em Link-Matrix} +\index{Link-Matrix}% +\begin{equation} H = \begin{pmatrix} -P(S'_1|S_1)&P(S'_1|S_2)&\dots &P(S'_1|S_N)\\ -P(S'_2|S_1)&P(S'_2|S_2)&\dots &P(S'_2|S_N)\\ +P(S'_1\mid S_1)&P(S'_1\mid S_2)&\dots &P(S'_1\mid S_N)\\ +P(S'_2\mid S_1)&P(S'_2\mid S_2)&\dots &P(S'_2\mid S_N)\\ \vdots &\vdots &\ddots&\vdots \\ -P(S'_N|S_1)&P(S'_N|S_2)&\dots &P(S'_N|S_N) +P(S'_N\mid S_1)&P(S'_N\mid S_2)&\dots &P(S'_N\mid S_N) \end{pmatrix}. -\] +\label{buch:google:eqn:linkmatrix} +\end{equation} Die Formel~\eqref{buch:google:eqn:totalewahrscheinlichkeit} wird dann zur Formel für das Produkt Matrix mal Vektor: \[ @@ -189,7 +201,7 @@ Formel für das Produkt Matrix mal Vektor: = \sum_{i=1}^N h_{ji} p_i = -\sum_{i=1}^N P(S'_j|S_i) P(S_i) +\sum_{i=1}^N P(S'_j\mid S_i) P(S_i) = p'_j \qquad\Rightarrow\qquad @@ -217,13 +229,26 @@ H = \end{equation} \qedhere \end{beispiel} - +Die Link-Matrix kann aus der Adjazenzmatrix des gerichteten Graphen +bestimmt werden. +Dazu ist zu beachten, dass jede Spalte durch die Anzahl der Einsen +in dieser Spalte zu teilen ist. +Ein Zeilenvektor, der die Zahl der Einsen enthält, entsteht durch +Multiplikation mit einem Zeilenvektor $U^t$ aus lauter Einsen. +Mit dem Hadamard-Produkt ist dann die Link-Matrix durch +\[ +H += +(U(U^tA(G))^{\odot(-1)})\odot A(G) +\] +gegeben, wobei $(U^tA(G))^{\odot(-1)}$ die Inverse bezüglich des +Hadamard-Produktes ist. % % Freier Wille % \subsection{``Freier Wille'' \label{buch:subsection:freier-wille}} -Das Modell in +Das in Abschnitt~\eqref{buch:subsection:wahrscheinlichkeitsinterpretation} beschriebene Modell geht unter anderem davon aus, dass der Benutzer ausschliesslich die Navigation entlang der Links verwendet. @@ -245,38 +270,39 @@ Wir bezeichnen das Ereignis, dass der Benutzer nicht die Link-Navigation verwendet mit $F$ für ``freier Wille'', obwohl es so etwas natürlich nicht gibt. Die Wahrscheinlichkeit, auf der Seite $S'_j$ zu landen, setzt sich jetzt -aus den zwei Fällen $F$ und $\overline{F}$ zusammen, für die erneut der +aus den zwei Fällen $F$ und $\smash{\overline{F}}$ zusammen, für die erneut der Satz von der totalen Wahrscheinlichkeit den Zusammenhang \[ P(S'_j) = -P(S'_j|\overline{F}) P(\overline{F}) +P(S'_j\mid \overline{F}) P(\overline{F}) + -P(S'_j|F) P(F) +P(S'_j\mid F) P(F) \] +liefert. Die Wahrscheinlichkeit $\alpha = P(F)$, mit der der Benutzer den -``freiene Willen'' bemüht, kann experimentell durch Studien ermittelt +``freien Willen'' bemüht, kann experimentell durch Studien ermittelt werden, die das Benutzerverhalten beobachten. -Die Wahrscheinlichkeit $P(S'_j|\overline{F})$ entsteht dadurch, dass +Die Wahrscheinlichkeit $P(S'_j\mid \overline{F})$ entsteht dadurch, dass der Benutzer der Linknavigation folgt, sie entspricht also der früher -berechnenten Wahrscheinlichkeit +berechneten Wahrscheinlichkeit \[ -P(S'_j|\overline{F}) = \sum_{i=1}^N P(S'_j|S_i) P(S_i). +P(S'_j\mid \overline{F}) = \sum_{i=1}^N P(S'_j\mid S_i) P(S_i). \] oder in Vektorform \[ -(P(S'_j|\overline{F}))_{j=1,\dots,n} +(P(S'_j\mid \overline{F}))_{j=1,\dots,n} = Hp. \] -Über die spontane Besuchswahrscheinlichkeit $P(S'_j|F)$ wissen wir +Über die spontane Besuchswahrscheinlichkeit $P(S'_j\mid F)$ wissen wir nichts. Eine erste Annahme könnte sein, dass jede Seite gleich wahrscheinlich -ist, dass also $P(S'_j|F)=1/N$. +ist, dass also $P(S'_j\mid F)=1/N$. Alternativ könnte man auch eine Wahrscheinlichkeitsverteilung -$q_j = P(S'_j|F)$ experimentell zu ermitteln versuchen. +$q_j = P(S'_j\mid F)$ experimentell zu ermitteln versuchen. Unter der Annahme, dass alle Seitenbesuche im Falle $F$ auf Grund eines Sucheresultats einer Suchmaschine erfolgen, könnte die Suchmaschine den Vektor $q$ aus ihrer eigenen Suchstatistik ermitteln. @@ -286,7 +312,7 @@ Das erweiterte Modell kann also durch P(S'_j) = \sum_{i=1}^N -\alpha P(S'_j|S_i) P(S_i) +\alpha P(S'_j\mid S_i) P(S_i) + (1-\alpha) q_j \qquad\Rightarrow\qquad @@ -309,7 +335,7 @@ Wir streben daher an, die Formel~\eqref{buch:google:eqn:composed} ebenfalls in die Form $p=Gp$ mit einer neuen Matrix $G$ zu bringen. Die Matrixform von -\label{buch:google:eqn:composed} +\eqref{buch:google:eqn:composed} zeigt, dass sich die gesuchte Matrix $G$ zusammensetzt aus dem Summanden $\alpha H$ und einem weiteren Summanden $A$ mit der Eigenschaft, dass $Ap = q$ für jeden beliebigen Wahrscheinlichkeitsvektor $p$. @@ -384,8 +410,10 @@ heisst die \index{Google-Matrix}% \end{definition} -Die Google-Matrix wurde von Sergei Brin und Larry Page -in dem Artikel \cite{BRIN1998107} als Basis der Suchmaschine +Die Google-Matrix wurde von Sergey Brin und Larry Page +\index{Brin, Sergey}% +\index{Page, Larry}% +in dem Artikel \cite{BRIN1998107} als Grundlage der Suchmaschine Google beschrieben. Sie war die Basis für den Erfolg von Google und wird dem Prinzip nach auch heute noch zur Rangierung der Suchresultate verwendet. @@ -406,6 +434,8 @@ ansehen. Die Kundenwebsite wird daher in den Suchresultaten weiter oben gezeigt. Das Problem rührt natürlich daher, dass alle Links als gleichermassen aussagekräftig betrachtet werden. +Solche Websites werden heutzutage von der Berechnung der Google-Matrix +ausgeschlossen. Die aktuell verwendete Variante der Google-Matrix ist natürlich ein Betriebsgeheimnis der Firma Google. @@ -417,7 +447,7 @@ Betriebsgeheimnis der Firma Google. \label{buch:subsection:wahrscheinlichkeitsverteilung}} Die Google-Matrix $G$ selbst interessiert weniger als die Wahrscheinlichkeitsverteilung $p$. -Ziel dieses Abschnittes, ist den Vektor $p$ zu berechnen. +Ziel dieses Abschnittes ist, den Vektor $p$ zu berechnen. \subsubsection{Stationäre Verteilung} Die Einträge $P(S_i)$ des Vektors $p$ geben die Wahrscheinlichkeit an, mit @@ -445,14 +475,17 @@ Gp = p. $p$ ist also ein Eigenvektor der Matrix $G$ zum Eigenwert $1$. Für ein sehr kleines Netzwerk wie im oben dargestellten Beispiel ist es -einfach, mit verbreiteten numerischen Algorithmen alle Eigenwerte und +einfach, mit gängigen numerischen Algorithmen alle Eigenwerte und Eigenvektoren zu finden. Benötigt wird allerdings nur der Eigenvektor zum Eigenwert $1$. \begin{beispiel} -Ein Eigenvektor zum Eigenwert $1$ der Matrix $G$, die aus der Matrix $H$ +Octave +\index{Octave} +findet den folgenden Eigenvektor zum Eigenwert $1$ der Matrix $G$, +die aus der Matrix $H$ von \eqref{buch:google:eqn:linkmatrixbeispiel} -und dem Vektor $q=\frac18u$ und $\alpha=0.9$ gebildet wurde, ist +und dem Vektor $q=\frac18U$ und $\alpha=0.9$ gebildet wurde: \[ p_0=\begin{pmatrix} 0.20100\\ @@ -491,10 +524,14 @@ erhält man die Wahrscheinlichkeitsverteilung $p$. \subsubsection{Potenzverfahren} -Die üblichen Algorithmen wie der Francis-Algorithmus zur Bestimmung -von Eigenwerten und Eigenvektoren ist für grosse Matrizen nicht praktikabel. +Die üblichen Algorithmen wie der von den meisten Softwarepaketen +verwendete Francis-Algorithmus \cite{francis:watkins_paper,buch:watkins} +\index{Francis-Algorithmus}% +zur Bestimmung von Eigenwerten +und Eigenvektoren ist für grosse Matrizen nicht praktikabel. Da aber $1$ der betragsgrösste Eigenwert ist, kann sehr oft ein zugehöriger Eigenvektor mit der nachfolgend beschriebenen {\em Potenzmethode} +\index{Potenzmethode}% gefunden werden. Sei $A$ eine $n\times n$-Matrix, der Einfachheit halber nehmen wir an, @@ -535,8 +572,8 @@ a_n\biggl(\frac{\lambda_n}{\lambda_1}\biggr)^k v_n. \] Da alle Brüche Betrag $<1$ haben, konvergiert die rechte Seite für $k\to\infty$ gegeben den ersten Summanden. -Durch wiederholte Anwendung von $A/\lambda_1$ auf einen (fast) belieibigen -Startvektor $v$ erhält man also eine Folge von Vektoren, die gegen eine +Durch wiederholte Anwendung von $A/\lambda_1$ auf einen (fast) beliebigen +Startvektor $v$ erhält man also eine Folge von Vektoren, die gegen einen Eigenvektor zum Eigenwert $\lambda_1$ konvergiert. Numerische Ungenauigkeiten können bewirken, dass die Iteration mit der diff --git a/buch/chapters/80-wahrscheinlichkeit/markov.tex b/buch/chapters/80-wahrscheinlichkeit/markov.tex index 0485714..1e30010 100644 --- a/buch/chapters/80-wahrscheinlichkeit/markov.tex +++ b/buch/chapters/80-wahrscheinlichkeit/markov.tex @@ -18,7 +18,10 @@ werden. % \subsection{Markov-Eigenschaft} % XXX Notation, Zustände, Übergangswahrscheinlichkeit -Ein stochastischer Prozess ist eine Familie von Zustandsvariablen +Ein stochastischer Prozess ist eine Familie von Zufallsvariablen +\index{stochastischer Prozess}% +\index{Prozess, stochastisch}% +\index{Zufallsvariable}% $X_t$ mit Werten in einer Menge $\mathcal{S}$ von Zuständen. Der Parameter $t$ wird üblicherweise als die Zeit interpretiert, er kann beliebige reelle Werte oder diskrete Werte annahmen, im letzten @@ -36,6 +39,7 @@ Zustands $s\in\mathcal{S}$ zu einem späteren Zeitpunkt $t_1>t_0$ zu studieren. Das Ereignis $\{X_t = x\}$ kann man sich als abhängig von der Vorgeschichte vorstellen. +\index{Vorgeschichte}% Die Vorgeschichte besteht dabei aus dem Eintreten gewisser Ereignisse \[ \{X_0=x_0\}, @@ -47,7 +51,7 @@ Die Vorgeschichte besteht dabei aus dem Eintreten gewisser Ereignisse zu früheren Zeiten $t_0<t_1<\dots<t_n<t$. Die bedingte Wahrscheinlichkeit \begin{equation} -P(X_t = x| +P(X_t = x \mid X_{t_n}=x_n\wedge X_{t_{n-1}}=x_{n-1}\wedge\dots\wedge X_{t_1}=x_1\wedge X_{t_0}=x_0) \label{buch:wahrscheinlichkeit:eqn:historybedingt} @@ -58,6 +62,7 @@ die Zustände $x_0,x_1,\dots,x_n$ durchlaufen hat. \subsubsection{Gedächtnislosigkeit} % XXX Gedächtnislösigkeit, Markov-Eigenschaft +\index{Markov-Eigenschaft}% In vielen Fällen ist nur der letzte durchlaufene Zustand wichtig. Die Zustände in den Zeitpunkten $t_0<\dots<t_{n-1}$ haben dann keinen Einfluss auf die Wahrscheinlichkeit. @@ -73,25 +78,26 @@ $x_0,\dots,x_n,x\in \mathcal{S}$ die Wahrscheinlichkeit~\eqref{buch:wahrscheinlichkeit:eqn:historybedingt} nicht von der Vorgeschichte abhängt, also \[ -P(X_t = x| +P(X_t = x\mid X_{t_n}=x_n\wedge X_{t_{n-1}}=x_{n-1}\wedge\dots\wedge X_{t_1}=x_1\wedge X_{t_0}=x_0) = -P(X_t = x| +P(X_t = x \mid X_{t_n}=x_n). \] \index{Markov-Eigenschaft} \end{definition} -Die Wahrscheinlichkeiten $P(X_t=x|X_s=y)$ mit $t>s$ bestimmen das +Die Wahrscheinlichkeiten $P(X_t=x\mid X_s=y)$ mit $t>s$ bestimmen das zeitliche Verhalten der Wahrscheinlichkeiten vollständig. Wir schreiben daher auch \[ p_{xy}(t, s) = -P(X_t = x|X_s=y) +P(X_t = x\mid X_s=y) \] für die sogenannte {\em transiente Übergangswahrscheinlichkeit}. +\index{transiente Übergangswahrscheinlichkeit}% Für eine endliche Menge von Zuständen, können die transienten Übergangswahrscheinlichkeiten auch als zeitabhängige quadratische Matrix $P(s,t)$ geschrieben werden, deren @@ -105,13 +111,14 @@ mit den Zuständen $x,y\in\mathcal{S}$ indiziert sind. \subsubsection{Die Chapman-Kolmogorov-Gleichung} % XXX Chapman-Kolmogorov-Gleichung +\index{Chapman-Kolmogorov-Gleichung}% Man beachte, dass in der Definition der Markov-Eigenschaft keine Voraussetzungen darüber gemacht werden, wie nahe am Zeitpunkt $t$ der letzte Zeitpunkt $t_n$ der Vorgeschichte liegt. Die transienten Übergangswahrscheinlichkeiten $p_{xy}(s,t)$ werden aber im allgemeinen davon abhängen, wie weit in der Vergangenheit der Zeitpunkt $s<t$ liegt. -Für eine näheren Zeitpunkt $\tau$ mit $s<\tau <t$ muss es daher +Für einen näheren Zeitpunkt $\tau$ mit $s<\tau <t$ muss es daher einen Zusammenhang zwischen den transienten Übergangswahrscheinlichkeiten $p_{xy}(s,\tau)$, $p_{xy}(\tau,t)$ und $p_{xy}(s,t)$ geben. @@ -187,16 +194,18 @@ Es ist üblich, für die Zeitpunkte ganze oder natürliche Zahlen zu verwenden. \begin{definition} -Eine diskrete Markov-Kette ist ein stochastischer Prozess +Eine {\em diskrete Markov-Kette} ist ein stochastischer Prozess $(X_t)_{t\in\mathbb{N}}$ mit Werten in $\mathcal{S}$, der die Markov-Eigenschaft \[ -P(X_{n+1}=x_{n+1}|X_n=x_n\wedge\dots X_0=x_0) +P(X_{n+1}=x_{n+1}\mid X_n=x_n\wedge\dots X_0=x_0) = -P(X_{n+1}=x_{n+1}|X_n=x_n) +P(X_{n+1}=x_{n+1}\mid X_n=x_n) \] hat. \end{definition} +\index{diskrete Markov-Kette}% +\index{Markov-Kette, diskret}% \begin{figure} \centering @@ -220,8 +229,9 @@ p_{11}(n+1,n) & \dots & p_{1s}(n+1,n)\\ p_{11}(n+1,n) & \dots & p_{1s}(n+1,n) \end{pmatrix}, \] -auch die $1$-Schritt Übergangswahrscheinlichkeit genannt, kann man jetzt +auch die $1$-Schritt-Übergangswahrscheinlichkeit genannt, kann man jetzt auch die Matrix der Überganswahrscheinlichkeiten für mehrere Schritte +\index{Ubergangswahrscheinlichkeit@Übergangswahrscheinlichkeit}% \[ T(n+m,n) = @@ -239,12 +249,12 @@ verwendet werden, wenn sie zwei Bedingungen erfüllt: \begin{enumerate} \item Die Einträge von $T$ müssen als Wahrscheinlichkeiten interpretiert werden können, sie müssen also alle zwischen $0$ und $1$ sein: -$0\le t_{ij}\le 1$ für $i,j\in\mathcal{S}$ +$0\le t_{i\!j}\le 1$ für $i,j\in\mathcal{S}$ \item Die Matrix muss alle möglichen Fälle erfassen. Dazu ist notwendig, dass sich die Wahrscheinlichkeiten aller Übergänge aus einem Zustand $j$ zu $1$ summieren, also \[ -\sum_{i\in\mathcal{S}} p_{ij} = 1. +\sum_{i\in\mathcal{S}} p_{i\!j} = 1. \] Die Summe der Elemente einer Spalte \end{enumerate} @@ -252,6 +262,7 @@ Die Summe der Elemente einer Spalte \begin{beispiel} Die Permutationsmatrix einer Permutation $\sigma\in S_n$ (Abschnitt~\label{buch:section:permutationsmatrizen}) +\index{Permutationsmatrix}% ist eine Matrix mit Einträgen $0$ und $1$, so dass die erste Bedingung erfüllt ist. In jeder Zeile oder Spalte kommt genau eine $1$ vor, so dass auch die @@ -269,8 +280,8 @@ p_i(n) = P(X_i=n) \] -geschrieben, die auch in einem Vektor $p(n)$ zusammengefasst -werden können. +geschrieben, die auch in einem Vektor $p(n)$ mit den Komponten +$p_i(n)$ zusammengefasst werden können. Die Matrix der Übergangswahrscheinlichkeiten erlaubt, die Verteilung $p(n+1)$ aus der Verteilung $p(n)$ zu berechnen. Nach dem Satz von der totalen Wahrscheinlichkeit ist nämlich @@ -278,9 +289,9 @@ Nach dem Satz von der totalen Wahrscheinlichkeit ist nämlich P(X_{n+1}=x) = \sum_{y\in\mathcal{S}} -P(X_{n+1}=x|X_n=y) P(X_n=y) +P(X_{n+1}=x\mid X_n=y) P(X_n=y) \qquad\text{oder}\qquad -p^{(n+1)} = T(n+1,n) p^{(n)} +p(n+1) = T(n+1,n) p(n) \] in Matrixform. Die Zeitentwicklung kann also durch Multiplikation mit der Übergangsmatrix @@ -288,6 +299,7 @@ berechnet werden. \subsubsection{Zeitunabhängige Übergangswahrscheinlichkeiten} % XXX Übergangswahrscheinlichkeit +\index{zeitunabhängige Übergangswahrscheinlichkeiten} Besonderes einfach wird die Situation, wenn die Übergangsmatrix $T(n+1,n)$ nicht von der Zeit abhängt. In diesem Fall ist $T(n+1,n) = T$ für alle $n$. @@ -311,32 +323,41 @@ homogene Markov-Kette mit Übergangsmatrix $T$, wenn $Tp=p$. \end{definition} Eine stationäre Verteilung ist offenbar ein Eigenvektor der Matrix -$T$ zum Eigenwert $1$. +$T$ zum Eigenwert $1$. Gefunden werden kann er als Lösung des Gleichungssystems $Tp=p$. -Dazu muss die Matrix $T-E$ singulär sein. -Die Summe einer Spalte von $T$ ist aber immer ein, da $E$ in jeder Spalte +Dazu muss aber die Matrix $T-I$ singulär sein, wie man wie folgt +einsehen kann. +Die Summe einer Spalte von $T$ ist aber immer $1$, da sich die +Wahrscheinlichkeiten zu $1$ summieren müssen. +Da die Einheitsmatrix $I$ in jeder Spalte genau eine $1$ enthält, ist die Summe der Einträge einer Spalte von -$T-E$ folglich $0$. -Die Summe aller Zeilen von $T-E$ ist also $0$, die Matrix $T-E$ +$I$ ebenfalls $1$. +Die Summe einer Spalte von $T-I$ ist folglich $0$. +Die Summe aller Zeilen von $T-I$ ist also $0$, die Matrix $T-I$ ist singulär. -Dies garantiert aber noch nicht, dass alle Einträge in diesem -Eigenvektor auch tatsächlich nichtnegativ sind. + +Dass $T-I$ singulär ist, garantiert aber noch nicht, +dass alle Einträge in einem zum Eigenwert $1$ +Eigenvektor auch tatsächlich nichtnegativ gewählt werden können. Die Perron-Frobienus-Theorie von +\index{Perron-Frobenius-Theorie}% Abschnitt~\ref{buch:section:positive-vektoren-und-matrizen} -beweist, dass sich immer ein Eigenvektor mit nichtnegativen -Einträgen finden lässt. +beweist, dass genau dies immer möglich ist. -Es ist aber nicht garantiert, dass eine stationäre Verteilung +Es ist nicht garantiert, dass eine stationäre Verteilung auch eindeutig bestimmt ist. Dieser Fall tritt immer ein, wenn die geometrische Vielfachheit des Eigenwerts $1$ grösser ist als $1$. In Abschnitt~\ref{buch:subsection:elementare-eigenschaften} werden Bedingungen an eine Matrix $T$ untersucht, die garantieren, -dass der Eigenraum zum Eigenvektor $1$ einedeutig bestimmt ist. +dass der Eigenraum zum Eigenvektor $1$ eindimensional ist. \begin{beispiel} -Als Beispiel dafür betrachten wir eine Permutation $\sigma\in S_n$ -und die zugehörige Permutationsmatrix $P$, +Als Beispiel dafür, dass der Eigenraum $\mathcal{E}_1(T)$ +mehrdimensional sein kann, betrachten wir eine Permutation $\sigma\in S_n$ +\index{Permutation}% +und die zugehörige Permutationsmatrix $P_\sigma$, +\index{Permutationsmatrix}% wie sie in Abschnitt~\label{buch:section:permutationsmatrizen} beschrieben worden ist. Wir verwenden die @@ -365,7 +386,8 @@ setzt. Die Konstruktion stellt sicher, dass sich die Komponenten zu $1$ summieren. Wir können aus dem Beispiel auch ableiten, dass die geometrische -Vielfachheit des Eigenvektors $1$ mindestens so gross ist wie die +Vielfachheit des Eigenwerts $1$ einer Permutationsmatrix $P_\sigma$ +mindestens so gross ist wie die Anzahl der Zyklen der Permutation $\sigma$. \end{beispiel} @@ -377,8 +399,9 @@ Die Zyklen können daher unabhängig voneinander studiert werden. Diese Idee kann auf allgemeine Markov-Ketten verallgemeinert werden. \begin{definition} -Zwei Zustände $i,j\in\mathcal{S}$ kommunizieren, wenn die -Übergangswahrscheinlichkeiten $T_{ij}(n) \ne 0$ und $T_{ij}(n)\ne 0$ sind +Zwei Zustände $i,j\in\mathcal{S}$ {\em kommunizieren}, wenn die +\index{kommunizieren}% +Übergangswahrscheinlichkeiten $T_{i\!j}(n) \ne 0$ und $T_{i\!j}(n)\ne 0$ sind für $n$ gross genug. \end{definition} @@ -407,12 +430,14 @@ Solche Markov-Ketten können unabhängig voneinander studiert werden. Die Bedingung der Irreduzibilität ist gleichbedeutend damit, dass für genügend grosses $n$ alle Matrixelemente von $T^n$ positiv sind. -Solche Matrizen nennt man positiv, +Solche Matrizen nennt man {\em positiv}, +\index{positive Matrix}% in Abschnitt~\ref{buch:section:positive-vektoren-und-matrizen} wird gezeigt, dass positive Matrizen immer eine eindeutige stationäre Verteilung haben. In Abbildung~\ref{buch:wahrscheinlichkeit:fig:markovzerfall} ist eine reduzible Markov-Kette dargestellt, die Zustandsmenge +\index{reduzible Markov-Kette}% zerfällt in zwei Teilmengen von Zuständen, die nicht miteinander kommunizieren. Ein irreduzible Markov-Kette liegt vor, wenn sich ähnlich wie @@ -420,7 +445,7 @@ in Abbildung~\ref{buch:wahrscheinlichkeit:fig:diskretemarkovkette} jeder Zustand von jedem anderen aus erreichen lässt. Wenn sich der Vektorraum $\mathbb{R}^n$ in zwei unter $T$ invariante -Unterräme zerlegen lässt, dann hat nach Wahl von Basen in den Unterräumen +Unterräume zerlegen lässt, dann hat nach Wahl von Basen in den Unterräumen die Matrix $T$ die Form \[ \left( @@ -483,7 +508,7 @@ Die stationären Verteilungen \operatorname{Stat}(T) = \{ -p\in\mathbb R_+^n\;|\; \text{$Tp=p $ und $\|p\|_1=1$} +p\in\mathbb R_+^n \mid \text{$Tp=p $ und $\|p\|_1=1$} \} \] bilden was man eine konvexe Menge nennt. @@ -495,7 +520,7 @@ Jede Verteilung auf der ``Verbindungsstrecke'' zwischen den beiden Verteilungen ist auch wieder stationär. \begin{definition} -Eine {\em konvexe Kombination} von Vektoren $v_1,\dots,v_k\in\mathbb{R^n}$ +Eine {\em konvexe Kombination} von Vektoren $v_1,\dots,v_k\in\mathbb{R}^n$ ist ein Vektor der Form \[ v=t_1v_1+\dots + t_kv_k @@ -512,7 +537,8 @@ wieder in $M$ ist. Die konvexen Kombinationen der Vektoren sind Linearkombination mit nichtnegativen Koeffizienten. Sie bilden im Allgemeinen -einen $(k-1)$-Simplex in $\mathbb{R}^n$. +einen $(k-1)$-Simplex in $\mathbb{R}^n$ (siehe auch +Abbildung~\ref{buch:wahrscheinlichkeit:fig:konvex}). Für zwei Punkte $x$ und $y$ bilden die konvexen Kombination $tx+(1-t)y$ für $t\in[0,1]$ die Verbindungsstrecke der beiden Vektoren. @@ -527,7 +553,7 @@ ihre Verbindungsstrecke enthält Im Beispiel der Google-Matrix wurde ein iterativer Algorithmus zur Berechnung des Pagerank verwendet. Es stellt sich daher die Frage, ob diese Methode für andere homogene -Markov-Ketten auch funkioniert. +Markov-Ketten auch funktioniert. Man beginnt also mit einer beliebigen Verteilung $p(0)$ und wendet die Übergangsmatrix $T$ wiederholt an. Es entsteht somit eine Folge $p(n) = T^np(0)$. @@ -546,8 +572,8 @@ Verteilung. Für eine stationäre Verteilung $p(0)$ ist die Folge $p(n)$ eine konstante Folge, sie konvergiert also gegen $p(0)$. Stationäre Verteilungen sind also automatisch Grenzverteilungen. -Falls der Raum der stationären Verteilungen mehrdimensional sind, -dann ist auch die Grenzverteilung nicht eindeutig bestimmt, selbst +Falls der Raum der stationären Verteilungen mehrdimensional ist, +braucht die Grenzverteilung nicht eindeutig bestimmt zu sein, selbst wenn sie existiert. Aber nicht einmal die Existenz einer Grenzverteilung ist garantiert, wie das folgende Beispiel zeigt. @@ -578,6 +604,8 @@ p(2)&=p(5)=p(8)=\dots =\begin{pmatrix}p_3(0)\\p_1(0)\\p_2(0)\end{pmatrix}. \end{align*} Die Folge $p(n)$ kann also nur dann konvergieren, wenn die drei Komponenten gleich sind. +Insbesondere gibt es keine Grenzverteilung, wenn sie nicht alle +gleich sind. \end{beispiel} \subsubsection{Erwartungswert und Varianz} @@ -588,11 +616,11 @@ zu berechnen. Dazu muss jedem Zustand ein Zahlenwert zugeordnet werden. Sei also \( -g: \mathcal{S}\to R +g: \mathcal{S}\to \mathbb{R} \) eine Funktion, die einem Zustand eine reelle Zahl zuordnet. Aus der Zufallsvariable $X_n$ des Zustands zur Zeit $n$ wird daraus -die Zufallsvariable $Y_n=g(X_n)$ des Wertes zur Zeit $n$. +die reellwertige Zufallsvariable $Y_n=g(X_n)$ des Wertes zur Zeit $n$. Die Abbildung $g$ kann auch als Vektor mit der Komponenten $g_i$ für $i\in\mathcal{S}$ betrachtet werden, wir verwenden für diesen Vektor wieder die Schreibweise $g$. @@ -634,7 +662,7 @@ definieren. In Abschnitt~\ref{buch:section:paradoxon-von-parrondo} wird ein Spiel vorgestellt, in dem der Gewinn davon abhängt, welcher Übergang stattfindet, nicht welcher Zustand erreicht wird. -Es git daher eine Matrix $G$ von Gewinnen, der Eintrag $g_{ij}$ ist +Es git daher eine Matrix $G$ von Gewinnen, der Eintrag $g_{i\!j}$ ist der Gewinn, der bei einem Übergang von Zustand $j$ in den Zustand $i$ ausgezahlt wird. Mit dieser Matrix lassen sich jetzt viele verschiedene Fragen beantworten: @@ -642,7 +670,7 @@ Mit dieser Matrix lassen sich jetzt viele verschiedene Fragen beantworten: \begin{frage} \label{buch:wahrscheinlichkeit:frage1} Mit welchem Gewinn kann man in Runde $n$ des Spiels rechnen, -wenn $p(n-1)$ die Verteilung zur Zeit $n-1$ ist? +wenn die Verteilung zur Zeit $n-1$ durch $p(n-1)$ gegeben ist? \end{frage} Der Erwartungswert ist @@ -664,15 +692,15 @@ einer Spielrunde im Zustand $i$ befindet? \end{frage} Dies ist der Spezialfall der Frage~\ref{buch:wahrscheinlichkeit:frage1} -für die Verteilung $p_j(n-1) = \delta_{ij}$. +für die Verteilung $p_j(n-1) = \delta_{i\!j}$. Der Erwartungswert ist die Summe der Spalte $j$ der Matrix $G\odot T$. Man kann das Produkt $U^t(G\odot T)$ also auch als eine Zeilenvektor von Gewinnerwartungen unter der Vorbedingung $X_{n-1}=j$ betrachten. \[ \begin{pmatrix} -E(Y|X_{n-1}=1) +E(Y\mid X_{n-1}=1) &\dots& -E(Y|X_{n-1}=n) +E(Y\mid X_{n-1}=n) \end{pmatrix} = U^t (G\odot T). @@ -681,6 +709,9 @@ Indem man $G$ durch $G^{\odot k}$ ersetzt, kann man beliebige höhere Momente berechnen. \subsection{Absorbierende Zustände} +In diesem Abschnitt gehen wir immer von einer irreduziblen Markov-Kette +aus. + % XXX Definition Eine Grenzverteilung beschreibt die relative Häufigkeit, mit der der Prozess in den verschiedenen Zuständen vorbeikommt. @@ -710,13 +741,13 @@ sie für alle zukünftigen Zustände in diesem Zustand. Eine Markov-Kette kann mehrere absorbierende Zustände haben, wie in Abbildung~\ref{buch:wahrscheinlichkeit:fig:abs} dargestellt. -Indem man die absorbierenden Zustände zuerst auflistet, bekommt die -Übergangsmatrix die Form +Indem man die absorbierenden Zustände zuerst auflistet, gefolgt von +den transienten Zustädnen, bekommt die Übergangsmatrix die Form \[ T= \left( \begin{array}{c|c} -E&R\\ +I&R\\ \hline 0&Q \end{array} @@ -732,7 +763,7 @@ T^2 = \left( \begin{array}{c|c} -E&R+RQ \\ +I&R+RQ \\ \hline 0&Q^2 \end{array} @@ -742,7 +773,7 @@ T^3 = \left( \begin{array}{c|c} -E&R+RQ+RQ^2 \\ +I&R+RQ+RQ^2 \\ \hline 0&Q^3 \end{array} @@ -754,18 +785,19 @@ T^k = \left( \begin{array}{c|c} -E&\displaystyle R\sum_{l=0}^{k-1} Q^l \\ +I&\displaystyle R\sum_{l=0}^{k-1} Q^l \\ \hline 0&Q^k \end{array} \right). \] -Da man früher oder später in einem absorbierenden Zustand landet, -muss $\lim_{k\to\infty} Q^k=0$ sein. +Wegen der angenommenen Irreduzibilität wird man +früher oder später in einem absorbierenden Zustand landet, +daher muss $\lim_{k\to\infty} Q^k=0$ sein. Die Summe in der rechten oberen Teilmatrix kann man als geometrische Reihe summieren, man erhält die Matrix \[ -\sum_{l=0}^{k-1} Q^l = (E-Q)^{-1}(E-Q^k), +\sum_{l=0}^{k-1} Q^l = (I-Q)^{-1}(I-Q^k), \] die für $k\to\infty$ gegen \[ @@ -773,7 +805,7 @@ N = \lim_{k\to\infty} \sum_{l=0}^{k-1} Q^l = -(E-Q)^{-1} +(I-Q)^{-1} \] konvergiert. Die Matrix $N$ heisst die {\em Fundamentalmatrix} der absorbierenden @@ -784,12 +816,13 @@ Markov-Kette. % XXX Absorptionszeit Wie lange dauert es im Mittel, bis der Prozess in einem Absorptionszustand $i$ stecken bleibt? +\index{Absorbtionszeit}% Die Fundamentalmatrix $N$ der Markov-Kette beantwortet diese Frage. -Wenn der Prozess genau im Schritt $k$ zum ersten Mal Zustand $i$ +Wenn der Prozess genau im Schritt $k$ zum ersten Mal im Zustand $i$ ankommt, dann ist $E(k)$ die mittlere Wartezeit. Der Prozess verbringt also zunächst $k-1$ Schritte in transienten -Zuständen, bevor er in einen absorbierenden Zustand wechselt. +Zuständen, bevor er in einen absorbierenden Zustand $i$ wechselt. Wir brauchen die Wahrscheinlichkeit für einen Entwicklung des Zustandes ausgehend vom Zustand $j$, die nach $k-1$ Schritten im Zustand $l$ @@ -808,7 +841,7 @@ innerhalb der Menge der Pfade, die auch tatsächlich absorbiert werden, das ist die bedingte Wahrscheinlichkeit \begin{equation} \begin{aligned} -P(X_k = i\wedge X_{k-1} = l \wedge X_0=j|X_k=i) +P(X_k = i\wedge X_{k-1} = l \wedge X_0=j\mid X_k=i) &= \frac{ P(X_k = i\wedge X_{k-1} = l \wedge X_0=j) @@ -833,25 +866,25 @@ E(k) &= \sum_{k=0}^\infty k( -q^{(k)}_{lj} +q^{(k)}_{l\!j} - -q^{(k-1)}_{lj} +q^{(k-1)}_{l\!j} ) \notag \\ &= \dots + -(k+1)( -q^{(k)}_{lj} +k( +q^{(k-1)}_{l\!j} - -q^{(k+1)}_{lj} +q^{(k)}_{l\!j} ) + -k( -q^{(k-1)}_{lj} +(k+1)( +q^{(k)}_{l\!j} - -q^{(k)}_{lj} +q^{(k+1)}_{l\!j} ) + \dots @@ -860,23 +893,44 @@ q^{(k)}_{lj} &= \dots + -q^{(k-1)}_{lj} +k +q^{(k-1)}_{l\!j} +\underbrace{ +\mathstrut +- +q^{(k)}_{l\!j} ++ +(k+1) +q^{(k)}_{l\!j} }_{\displaystyle q^{(k)}_{l\!j}} +\mathstrut +- +(k+1) +q^{(k+1)}_{l\!j} ++ +\dots +\\ +&= +\dots ++ +q^{(k)}_{l\!j} + \dots = -\sum_{k} q^{(k)}_{lj}. +\sum_{k} q^{(k)}_{l\!j}. \notag \end{align} In zwei benachbarten Termen in \eqref{buch:wahrscheinlichkeit:eqn:telescope} -heben sich die Summanden $kq^{(k)}_{lj}$ weg, man spricht von +heben sich die Summanden $kq^{(k)}_{l\!j}$ weg, man spricht von einer teleskopischen Reihe. +\index{teleskopische Reihe}% Die verbleibenden Terme sind genau die Matrixelemente der Fundamentalmatrix $N$. Die Fundamentalmatrix enthält also im Eintrag $(l,j)$ die Wartezeit bis zur Absorption über den Zustand $l$. \subsubsection{Wartezeit} % XXX Mittlere Zeit bis zu einem bestimmten Zustand +\index{Wartezeit}% Die mittlere Wartezeit bis zum Erreichen eines Zustands kann mit der Theorie zur Berechnung der Absorptionszeit berechnet werden. Dazu modifiziert man den Prozess dahingehend, dass der Zielzustand diff --git a/buch/chapters/80-wahrscheinlichkeit/parrondo.tex b/buch/chapters/80-wahrscheinlichkeit/parrondo.tex index 50e7fda..94b39fc 100644 --- a/buch/chapters/80-wahrscheinlichkeit/parrondo.tex +++ b/buch/chapters/80-wahrscheinlichkeit/parrondo.tex @@ -32,7 +32,7 @@ E(X) = 1\cdot P(X=1) + (-1)\cdot P(X=-1) = -\frac12+e + (-1)\biggl(\frac12-e\biggr) +\frac12+e + (-1)(\frac12-e) = 2e. \) @@ -41,6 +41,7 @@ Die Gewinnerwartung ist also genau dann negativ, wenn $e<0$ ist. \subsubsection{Das Spiel $B$} Das zweite Spiel $B$ ist etwas komplizierter, da der Spielablauf vom aktuellen Kapital $K$ des Spielers abhängt. +\index{Kapital}% Wieder gewinnt oder verliert der Spieler eine Einheit, die Gewinnwahrscheinlichkeit hängt aber vom Dreierrest des Kapitals ab. Sei $Y$ die Zufallsvariable, die den Gewinn beschreibt. @@ -49,9 +50,9 @@ andernfalls ist sie $\frac34$. Formell ist \begin{equation} \begin{aligned} -P(Y=1|\text{$K$ durch $3$ teilbar}) &= \frac{1}{10} +P(Y=1\mid \text{$K$ durch $3$ teilbar}) &= \frac{1}{10} \\ -P(Y=1|\text{$K$ nicht durch $3$ teilbar}) &= \frac{3}{4} +P(Y=1\mid \text{$K$ nicht durch $3$ teilbar}) &= \frac{3}{4} \end{aligned} \label{buch:wahrscheinlichkeit:eqn:Bwahrscheinlichkeiten} \end{equation} @@ -74,7 +75,7 @@ statt, der Eintrag $b_{ij}$ ist die Wahrscheinlichkeit \[ b_{ij} = -P(K\equiv i|K\equiv j), +P(K\equiv i\mid K\equiv j), \] dass ein Übergang vom Zustand $j$ in den Zustand $i$ stattfindet. Die Matrix ist @@ -95,11 +96,11 @@ Mit den Wahrscheinlichkeiten von findet man die Gewinnerwartung \begin{equation} \begin{aligned} -E(Y| \text{$K$ durch $3$ teilbar}) +E(Y\mid \text{$K$ durch $3$ teilbar}) &= -1\cdot P(Y=1|K\equiv 0\mod 3) +1\cdot P(Y=1\mid K\equiv 0\mod 3) + -(-1)\cdot P(Y=-1|K\equiv 0\mod 3) +(-1)\cdot P(Y=-1\mid K\equiv 0\mod 3) \\ &= \frac1{10} @@ -108,11 +109,11 @@ E(Y| \text{$K$ durch $3$ teilbar}) = -\frac{8}{10} \\ -E(Y| \text{$K$ nicht durch $3$ teilbar}) +E(Y\mid \text{$K$ nicht durch $3$ teilbar}) &= -1\cdot P(Y=1|K\not\equiv 0\mod 3) +1\cdot P(Y=1\mid K\not\equiv 0\mod 3) + -(-1)\cdot P(Y=-1|K\not\equiv 0\mod 3) +(-1)\cdot P(Y=-1\mid K\not\equiv 0\mod 3) \\ &= \frac34-\frac14 @@ -131,9 +132,9 @@ Die Gewinnerwartung in diesem Fall ist dann \begin{align} E(Y) &= -E(Y|\text{$K$ durch $3$ teilbar}) \cdot \frac13 +E(Y\mid \text{$K$ durch $3$ teilbar}) \cdot \frac13 + -E(Y|\text{$K$ nicht durch $3$ teilbar}) \cdot \frac23 +E(Y\mid\text{$K$ nicht durch $3$ teilbar}) \cdot \frac23 \notag \\ &= @@ -164,13 +165,13 @@ G=\begin{pmatrix} \end{pmatrix} \] gibt die Gewinne an, die bei einem Übergang anfallen. -Die Matrixelemente $g_{ij}b_{ij}$ des Hadamard-Produktes -$G\odot B$ -von $G$ mit $B$ enthält in den Spalten die Gewinnerwartungen +Die Matrix mit den Matrixelementen $g_{ij}b_{ij}$ ist das Hadamard-Produktes +$G\odot B$ von $G$ mit $B$. +Sie enthält in den Spalten die Gewinnerwartungen für die einzelnen Übergänge aus einem Zustand. Die Summe der Elemente der Spalte $j$ enthält die Gewinnerwartung \[ -E(Y|K\equiv j) +E(Y\mid K\equiv j) = \sum_{i=0}^2 g_{ij}b_{ij} \] @@ -181,9 +182,9 @@ $U^t=\begin{pmatrix}1&1&1\end{pmatrix}$ entsteht: \[ \begin{pmatrix} -E(Y|K\equiv 0)& -E(Y|K\equiv 1)& -E(Y|K\equiv 2) +E(Y\mid K\equiv 0)& +E(Y\mid K\equiv 1)& +E(Y\mid K\equiv 2) \end{pmatrix} = U^t @@ -194,7 +195,7 @@ Die Gewinnerwartung ist dann das Produkt E(Y) = \sum_{i=0}^2 -E(Y|K\equiv i) p_i +E(Y\mid K\equiv i) p_i = U^t (G\odot B)p. @@ -247,7 +248,7 @@ Das Spiel kennt die Dreierreste als die drei für das Spiel ausschlaggebenden Zuständen. Das Zustandsdiagramm~\ref{buch:wahrscheinlichkeit:fig:spielB} zeigt die möglichen Übergänge und ihre Wahrscheinlichkeiten, die zugehörige -Matrix ist +Übergangsmatrix ist \[ B = @@ -255,7 +256,7 @@ B 0 &\frac14 &\frac34\\ \frac1{10} &0 &\frac14\\ \frac9{10} &\frac34 &0 -\end{pmatrix} +\end{pmatrix}. \] Die Matrix $B$ ist nicht negativ und man kann nachrechnen, dass $B^2>0$ ist. Damit ist die Perron-Frobenius-Theorie von @@ -263,6 +264,7 @@ Abschnitt~\ref{buch:section:positive-vektoren-und-matrizen} anwendbar. Ein Eigenvektor zum Eigenwert $1$ kann mit Hilfe des Gauss-Algorithmus +\index{Gauss-Algorithmus}% gefunden werden: \begin{align*} \begin{tabular}{|>{$}c<{$}>{$}c<{$}>{$}c<{$}|} @@ -301,7 +303,7 @@ Daraus liest man einen möglichen Lösungsvektor mit den Komponenten $5$, $2$ und $6$ ab. Wir suchen aber einen Eigenvektor, der als Wahrscheinlichkeitsverteilung dienen kann. -Dazu müssen sich die Komponente zu $1$ summieren, was man durch normieren +Dazu müssen sich die Komponenten zu $1$ summieren, was man durch Normieren in der $l^1$-Norm erreichen kann: \begin{equation} p @@ -344,11 +346,13 @@ nach \begin{align*} P(Y=+1) &= -P(Y=+1|K\equiv 0) \cdot P(K\equiv 0) +P(Y=+1\mid K\equiv 0) \cdot P(K\equiv 0) + -P(Y=+1|K\equiv 1) \cdot P(K\equiv 1) +P(Y=+1\mid K\equiv 1) \cdot P(K\equiv 1) +\\ +&\qquad + -P(Y=+1|K\equiv 2) \cdot P(K\equiv 2) +P(Y=+1\mid K\equiv 2) \cdot P(K\equiv 2) \\ &= \frac{1}{10}\cdot\frac{5}{13} @@ -368,11 +372,13 @@ P(Y=+1|K\equiv 2) \cdot P(K\equiv 2) \\ P(Y=-1) &= -P(Y=-1|K\equiv 0) \cdot P(K\equiv 0) +P(Y=-1\mid K\equiv 0) \cdot P(K\equiv 0) + -P(Y=-1|K\equiv 1) \cdot P(K\equiv 1) +P(Y=-1\mid K\equiv 1) \cdot P(K\equiv 1) +\\ +&\qquad + -P(Y=-1|K\equiv 2) \cdot P(K\equiv 2) +P(Y=-1\mid K\equiv 2) \cdot P(K\equiv 2) \\ &= \frac{9}{10}\cdot\frac{5}{13} @@ -479,9 +485,9 @@ G\odot F = \begin{pmatrix} Nach der früher dafür gefundenen Formel ist \begin{align*} \begin{pmatrix} -E(Y|K\equiv 0)& -E(Y|K\equiv 1)& -E(Y|K\equiv 2) +E(Y\mid K\equiv 0)& +E(Y\mid K\equiv 1)& +E(Y\mid K\equiv 2) \end{pmatrix} &= U^t (G\odot \tilde{B}) @@ -710,10 +716,10 @@ A=\begin{pmatrix} \subsubsection{Das Spiel $C$} In jeder Durchführung des Spiels wird mit einem Münzwurf entschieden, ob Spiel $A$ oder Spiel $B$ gespielt werden soll. -Mit je Wahrscheinlichkeit $\frac12$ werden also die Übergansmatrizen +Mit Wahrscheinlichkeit je $\frac12$ werden also die Übergansmatrizen $A$ oder $B$ verwendet: \[ -P(K\equiv i|K\equiv j) +P(K\equiv i\mid K\equiv j) = A\cdot P(\text{Münzwurf Kopf}) + diff --git a/buch/chapters/80-wahrscheinlichkeit/positiv.tex b/buch/chapters/80-wahrscheinlichkeit/positiv.tex index 9f8f38f..159d6d3 100644 --- a/buch/chapters/80-wahrscheinlichkeit/positiv.tex +++ b/buch/chapters/80-wahrscheinlichkeit/positiv.tex @@ -7,15 +7,18 @@ \label{buch:section:positive-vektoren-und-matrizen}} \rhead{Positive Vektoren und Matrizen} Die Google-Matrix und die Matrizen, die wir in Markov-Ketten angetroffen +\index{Google-Matrix}% haben, zeichnen sich dadurch aus, dass alle ihre Einträge positiv oder mindestens nicht negativ sind. Die Perron-Frobenius-Theorie, die in diesem Abschnitt entwickelt +\index{Perron-Frobenius-Theorie}% werden soll, zeigt, dass Positivität einer Matrix nützliche Konsequenzen für Eigenwerte und Eigenvektoren hat. -Das wichtigste Resultat ist die Tatsache, dass postive Matrizen immer +Das wichtigste Resultat ist die Tatsache, dass positive Matrizen immer einen einzigen einfachen Eigenwert mit Betrag $\varrho(A)$ haben, -was zum Beispiel die Konvergenz des Pagerank-Algorithmus garantiert. -Dies wird im Satz von Perron-Frobenius in +was zum Beispiel die Konvergenz des PageRank-Algorithmus garantiert. +Dies wird im Satz~\ref{buch:wahrscheinlichkeit:satz:perron-frobenius} +von Perron-Frobenius in Abschnitt~\ref{buch:subsection:der-satz-von-perron-frobenius} erklärt. @@ -42,6 +45,8 @@ seine Komponenten nicht negativ sind: $v_i\ge 0\forall i$. Geometrisch kann man sich die Menge der positven Vektoren in zwei Dimensionen als die Punkte des ersten Quadranten oder in drei Dimensionen als die +\index{Quadrant}% +\index{Oktant}% Vektoren im ersten Oktanten vorstellen. Aus der Positivität eines Vektors lässt sich jetzt eine Vergleichsrelation @@ -62,9 +67,9 @@ Die Definition funktionieren analog auch für Matrizen: \begin{definition} Eine Matrix $A\in M_{m\times n}(\mathbb{R})$ heisst {\em positiv}, -wenn alle ihre Einträge $a_{ij}$ positiv sind: $a_{ij}>0\forall i,j$. +wenn alle ihre Einträge $a_{i\!j}$ positiv sind: $a_{i\!j}>0\forall i,j$. Eine Matrix $A\in M_{m\times n}(\mathbb{R})$ heisst {\em nichtnegativ}, -wenn alle ihre Einträge $a_{ij}$ nichtnegativ sind: $a_{ij}\ge 0\forall i,j$. +wenn alle ihre Einträge $a_{i\!j}$ nichtnegativ sind: $a_{i\!j}\ge 0\forall i,j$. \index{positive Matrix}% \index{nichtnegative Matrix}% Man schreibt $A>B$ bzw.~$A\ge B$ wenn $A-B>0$ bzw.~$A-B\ge 0$. @@ -133,7 +138,7 @@ Die Zyklen-Zerlegung einer Permutationsmatrix zeigt, welche Unterräume von $\mathbb{R}^n$ die iterierten Bilder eines Standardbasisvektors aufspannen. Diese sind invariante Unterräume der Matrix. -Das im Beispiel illustrierte Phänomen findet dann nur in invarianten +Das im Beispiel illustrierte Phänomen findet nur in invarianten Unterräumen statt. \begin{beispiel} @@ -151,7 +156,7 @@ A=\begin{pmatrix} \end{equation} besteht aus zwei $3\times 3$-Blöcken. Die beiden Unterräume $V_1=\langle e_1,e_2,e_3\rangle$ -und $V_2=\langle e_4,e_5,e_6\rangle$ sind daher invariante +und $V_2=\langle e_4,e_5,e_6\rangle$ sind invariante Unterräume von $A$ und damit auch von $A^n$. Die Potenzen haben daher auch die gleich Blockstruktur. Insbesondere sind zwar die Blöcke von $A^n$ für $n>1$ positive @@ -161,6 +166,7 @@ Teilmatrizen, aber die Matrix $A^n$ ist für alle $n$ nicht positiv. \begin{definition} Eine nichtnegative Matrix mit der Eigenschaft, dass $A^n>0$ für ein genügend grosses $n$, heisst {\em primitiv}. +\index{primitive Matrix}% \end{definition} Die Matrix $A$ von \eqref{buch:wahrscheinlichkeit:eqn:diffusion} @@ -200,6 +206,7 @@ hinaus. \begin{satz}[Trenntrick] \label{buch:wahrscheinlichkeit:satz:trenntrick} +\index{Trenntrick}% Sind $u$ und $v$ nichtnegative Vektoren und $u>v$, dann gibt es eine positive Zahl $\varepsilon>0$ derart, dass $u\ge (1+\varepsilon)v$. @@ -214,7 +221,7 @@ Wir betrachten die Zahl = \max_{v_i\ne 0} \frac{u_i}{v_i}. \] -Wegen $u>v$ sind die Quotienten auf der rechten Seite alle $>0$. +Wegen $u>v$ sind die Quotienten auf der rechten Seite alle $>1$. Da nur endlich viele Quotienten miteinander verglichen werden, ist daher auch $\vartheta >1$. Es folgt $u\ge \vartheta v$. @@ -244,6 +251,7 @@ $Au>Av$ (siehe auch Satz~\ref{buch:wahrscheinlichkeit:satz:vergleichstrick}) \begin{satz}[Vergleichstrick] \label{buch:wahrscheinlichkeit:satz:vergleichstrick} +\index{Vergleichstrick}% Sei $A$ eine positive Matrix und seinen $u$ und $v$ Vektoren mit $u\ge v$ und $u\ne v$, dann ist $Au > Av$ (siehe auch Abbildung~\ref{buch:wahrscheinlichkeit:fig:vergleich}). @@ -251,8 +259,8 @@ mit $u\ge v$ und $u\ne v$, dann ist $Au > Av$ \begin{proof}[Beweis] Wir schreiben $d=u-v$, nach Voraussetzung ist $d\ne 0$. -Der Satz besagt dann, dass aus $d\ge 0$ folgt, dass $Ad>0$, dies -müssen wir beweisen. +Der Satz besagt dann, dass aus $d\ge 0$ folgt, dass $Ad>0$. +Dies müssen wir beweisen. Die Ungleichung $Ad>0$ besagt, dass alle Komponenten von $Ad$ positiv sind. @@ -261,16 +269,16 @@ Um dies nachzuweisen, berechnen wir (Ad)_i = \sum_{j=1}^n -a_{ij} +a_{i\!j} d_j. \label{buch:wahrscheinlichkeit:eqn:Adpositiv} \end{equation} -Alle Terme $a_{ij}>0$, weil $A$ positiv ist, und mindestens eine -der Komponenten $d_j>0$, weil $d\ne 0$. +Alle Koeffizienten $a_{i\!j}$ sind $>0$, weil $A$ positiv ist. +Mindestens eine der Komponenten $d_j$ ist $>0$, weil $d\ne 0$. Insbesondere sind alle Terme der Summe $\ge 0$, woraus wir bereits schliessen können, dass $(Ad)_i\ge 0$ sein muss. Die Komponente $d_j>0$ liefert einen positiven Beitrag -$a_{ij}d_j>0$ +$a_{i\!j}d_j>0$ zur Summe~\eqref{buch:wahrscheinlichkeit:eqn:Adpositiv}, also ist $(Ad)_i>0$. \end{proof} @@ -284,8 +292,8 @@ Ist $A$ eine positive Matrix und $u\ge 0$ mit $u\ne 0$, dann ist $Au>0$. \end{korollar} -Eine positive Matrix macht also aus nicht verschwindenden -und nicht negativen Vektoren positive Vektoren. +Eine positive Matrix macht also aus nicht verschwindenden, +nicht negativen Vektoren positive Vektoren. % % Die verallgemeinerte Dreiecksungleichung @@ -331,30 +339,38 @@ eines gemeinsamen Einheitsvektors $c$ sind: $u_i=|u_i|c$ \begin{proof}[Beweis] Die Aussage kann mit vollständiger Induktion bewiesen werden. -Die Induktionsverankerung ist der Fall $n=2$ gegeben durch die +Die Induktionsverankerung ist der Fall $n=2$, gegeben durch die gewöhnliche Dreiecksungleichung. Wir nehmen daher jetzt an, die Aussage sei für $n$ bereits bewiesen, -wir müssen sie dann für $n+1$ beweisen. +wir müssen sie für $n+1$ beweisen. Die Summe von $n+1$ Vektoren kann man $u=u_1+\dots+u_n$ und $v=u_{n+1}$ aufteilen. -Es gilt dann +Es gilt nach der gewöhnlichen Dreiecksungleichung, dass \[ |u+v| = |u_1+\dots+u_n+u_{n+1}| +\le +|u_1+\dots+u_n|+|u_{n+1}| \] -und +mit Gleichheit genau dann, wenn $u_1+\dots+u_n$ und $u_{n+1}$ +linear abhängig sind. +Nach Induktionsannahme gilt ausserdem \[ -|u_1+\dots+u_n| = |u_1|+\dots+|u_n|. +|u_1+\dots+u_n| \le |u_1|+\dots+|u_n| \] -Aus der Induktionsannahme folgt dann, dass die Vektoren $u_1,\dots,u_n$ +mit Gleichheit genau dann, wenn die Vektoren $u_1,\dots,u_n$ positive Vielfache eines Einheitsvektors $u$ sind, $u_i=|u_i|c$. Es ist dann \[ -u=u_1+\dots+u_n = \biggl(\sum_{i=1}^n |u_i|\biggr). +u=u_1+\dots+u_n += +\biggl(\sum_{i=1}^n |u_i|c\biggr) += +\biggl(\sum_{i=1}^n |u_i|\biggr)c. \] -Aus der gewöhnlichen Dreiecksungleichung, angewendet auf $u$ und $v$ +Da $|u+v|=|u|+|v|$ genau dann gilt, wenn $u$ und $v$ linear abhängig sind, folgt jetzt, dass $v$ ebenfalls ein nichtnegatives Vielfaches von $c$ ist. Damit ist der Induktionsschritt vollzogen. \end{proof} @@ -380,7 +396,7 @@ Die motiviert den nachstehenden geometrischen Beweis des Satzes. \begin{proof}[Beweis] Wer stellen uns die komplexen Zahlen $u_i$ als Vektoren in der -zweidimensionalen Gaussschen Ebene vor. +zweidimensionalen Gauss\-schen Ebene vor. Dann ist die Aussage nichts anderes als ein Spezialfall von Satz~\ref{buch:wahrscheinlichkeit:satz:verallgemeinerte-dreiecksungleichung} für den zweidimensionalen reellen Vektorraum $\mathbb{C}$. @@ -396,8 +412,8 @@ Wir sind an den Eigenwerten und Eigenvektoren einer positiven oder primitiven Matrix interessiert. Nach Definition des Spektralradius $\varrho(A)$ muss es einen Eigenvektor zu einem Eigenwert $\lambda$ mit Betrag $|\lambda|=\varrho(A)$ geben, -aber a priori wissen wir nicht, ob es einen reellen Eigenwert vom -Betrag $\varrho(A)$ gibt, und ob der Eigenvektor dazu reell ist. +aber a priori wissen wir nicht, ob es einen reellen Eigenvektor zum +Eigenwert $\varrho(A)$ gibt. \begin{figure} \centering @@ -415,14 +431,16 @@ Iteriert man dies (Abbildung~\ref{buch:wahrscheinlichkeit:figure:positiv}), wird die Bildmenge immer enger, bis sie nur ein sehr enger Kegel um die Richtung des Eigenvektors ist. Tatsächlich kann man aus dieser Idee auch einen topologischen -Beweis des untenstehenden Satzes von Perron-Frobenius konstruieren. +Beweis des untenstehenden Satzes von Perron-Frobenius konstruieren +(\cite{skript:pftopo} und +\cite{skript:hilbertmetric}). Er beruht darauf, dass eine Abbildung, die Distanzen verkleinert, einen Fixpunkt hat. Die Konstruktion einer geeigneten Metrik ist allerdings eher kompliziert, weshalb wir im Beweise der nachstehenden Aussagen den konventionellen Weg wählen. -Wir beginnen damit zu zeigen, dass für positive Matrizen $A$, +Wir beginnen damit zu zeigen, dass für positive Matrizen $A$ nichtnegative Eigenvektoren zu Eigenwerten $\lambda\ne 0$ automatisch positiv sind. Ausserdem müssen die zugehörigen Eigenwerte sogar positiv sein. @@ -444,6 +462,14 @@ alle Komponenten von $\lambda u$ positiv sein. Das ist nur möglich, wenn $\lambda > 0$. \end{proof} +Wenn $v$ ein Eigenvektor von $A$ ist, dann ist auch jedes Vielfache +davon ein Eigenvektor, insbesondere können einzelne Komponenten +des Vektors $v$ auch negativ sein. +Der folgende Satz zeigt aber, dass man der Vektor aus den Beträgen +von der Komponenten von $v$ ebenfalls ein Eigenvektor zum +gleichen Eigenwert ist. +Insbesondere gibt es immer einen nichtnegativen Eigenvektor. + \begin{satz} \label{buch:wahrscheinlichkeit:satz:positivereigenvektor} Sei $A$ eine positive Matrix und $v$ ein Eigenvektor von $A$ zu einem @@ -457,12 +483,12 @@ Es gilt natürlich auch, dass \[ (Au)_i = -\sum_{j=1}^n a_{ij}u_j +\sum_{j=1}^n a_{i\!j}u_j = -\sum_{j=1}^n |a_{ij}v_j| +\sum_{j=1}^n |a_{i\!j}v_j| \ge \biggl| -\sum_{j=1}^n a_{ij}v_j +\sum_{j=1}^n a_{i\!j}v_j \biggr| = |(Av)_i| @@ -490,22 +516,26 @@ können wir jetzt eine Zahl $\vartheta>1$ finden derart, dass A^2 u \ge \vartheta \varrho(A) Au \] ist. -Durch weitere Anwendung von $A$ findet man -\begin{align*} +Durch wiederholte Anwendung von $A$ findet man +\begin{align} A^3 u & \ge (\vartheta \varrho(A))^2 Au +\notag \\ &\phantom{0}\vdots +\notag \\ A^{k+1} u & \ge (\vartheta \varrho(A))^{k} Au -\end{align*} -Daraus kann man jetzt die Norm abschätzen: +\label{buch:pf:eqn:ak+1} +\end{align} +Aus $|A^{k+1}u| \le \|A^k\|\,|Ak|$ und +\eqref{buch:pf:eqn:ak+1} kann man jetzt die Norm von $A^k$ abschätzen: \[ \begin{aligned} -\| A^{k}\|\, |Au| +\| A^{k}\|\cdot |Au| &\ge -\| A^{k+1}u\| +| A^{k+1}u| \ge -(\vartheta\varrho(A))^{k} |Au| +(\vartheta\varrho(A))^{k}\, |Au| && \Rightarrow & @@ -518,8 +548,11 @@ Daraus kann man jetzt die Norm abschätzen: \lim_{k\to\infty} \|A^k\|^{\frac{1}{k}} &\ge \vartheta\varrho(A) \\ -&&&\Rightarrow& -\varrho(A)&\ge \vartheta\varrho(A) +&&&&\|\phantom{00}& +\\ +&&&%\Rightarrow& +& +\varrho(A)&\ge \vartheta\varrho(A). \end{aligned} \] Wegen $\vartheta>1$ ist dies aber gar nicht möglich. @@ -527,6 +560,10 @@ Dieser Widerspruch zeigt, dass $u=v$ sein muss, insbesondere ist $v$ ein nichtnegativer Eigenvektor. \end{proof} +Die Potenzmethode funktioniert nur, wenn kein anderer Eigenwert +den Betrag $\varrho(A)$ hat. +Der folgende Satz garantiert dies. + \begin{satz} Sei $A$ eine positive Matrix und $v$ ein Eigenvektor zu einem Eigenwert $\lambda$ mit Betrag $|\lambda|=\varrho(A)$. @@ -541,24 +578,22 @@ Aus der Eigenvektorgleichung für $u$ folgt \begin{equation} Au = \varrho(A) u \quad\Rightarrow\quad -\sum_{j=1}^n a_{ij}|v_j| = \varrho(A) |v_i|. +\sum_{j=1}^n a_{i\!j}|v_j| = \varrho(A) |v_i|. \label{buch:wahrscheinlichkeit:eqn:pev1} \end{equation} Anderseits ist $v$ ein Eigenvektor zum Eigenwert $\lambda$, also gilt \[ -\sum_{j=1}^n a_{ij}v_j = \lambda v_i. +\sum_{j=1}^n a_{i\!j}v_j = \lambda v_i. \] Der Betrag davon ist \begin{equation} \biggl| -\sum_{j=1}^n a_{ij}v_j +\sum_{j=1}^n a_{i\!j}v_j \biggr| = |\lambda v_i| = -\varrho(A) |v_i| -= -\varrho |v_i|. +\varrho(A) |v_i|. \label{buch:wahrscheinlichkeit:eqn:pev2} \end{equation} Die beiden Gleichungen @@ -566,28 +601,35 @@ Die beiden Gleichungen und \eqref{buch:wahrscheinlichkeit:eqn:pev2} zusammen ergeben die Gleichung -\[ +\begin{equation} \biggl| -\sum_{j=1}^n a_{ij}v_j +\sum_{j=1}^n a_{i\!j}v_j \biggr| = -\sum_{j=1}^n a_{ij}|v_j|. -\] +\sum_{j=1}^n a_{i\!j}|v_j|. +\label{buch:pf:eqn:gleich} +\end{equation} Nach der verallgemeinerten Dreiecksungleichung Satz~\ref{buch:subsection:verallgemeinerte-dreiecksungleichung} -folgt jetzt, dass es eine komplexe Zahl $c$ vom Betrag $1$ gibt derart, +folgt jetzt aus der Gleichheit in~\eqref{buch:pf:eqn:gleich}, +dass es eine komplexe Zahl $c$ vom Betrag $1$ gibt derart, dass $v_j = |v_j|c=u_jc$. -Insbesondere ist $v=cu$ und damit ist +Insbesondere ist $v=cu$. +Damit kann man jetzt $\lambda$ berechnen, es ist \[ \lambda v = Av = Acu = c Au = c\varrho(A) u = \varrho(A) v, \] woraus $\lambda=\varrho(A)$ folgt. \end{proof} +In Anwendungen wollen wir schliessen, dass die Grenzverteilung +eindeutig ist, dazu ist notwendig, dass der Eigenraum des +Eigenwertes $\varrho(A)$ eindimensional ist. + \begin{satz} \label{buch:wahrscheinlichkeit:satz:geometrischeinfach} -Der Eigenraum einer positiven Matrix $A$ zum Eigenwert $\varrho(A)$ ist -eindimensional. +Der Eigenraum $E_{\varrho(A)}(A)$ einer positiven Matrix $A$ +zum Eigenwert $\varrho(A)$ ist eindimensional. \end{satz} \begin{proof}[Beweis] @@ -613,7 +655,7 @@ A(u-cv) \] Der Vektor auf der rechten Seite hat mindestens eine verschwindende Komponente. -Der Vektor auf der linken Seite ist nach Vergleichstrick +Der Vektor auf der linken Seite ist nach dem Vergleichstrick Satz~\ref{buch:wahrscheinlichkeit:satz:vergleichstrick} \[ A(u-cv) > 0, @@ -623,14 +665,20 @@ Dieser Widerspruch zeigt, dass die Annahme, es gäbe einen von $u$ linear unabhängigen Eigenvektor zum Eigenwert $\varrho(A)$ nicht haltbar ist. \end{proof} +Satz~\ref{buch:wahrscheinlichkeit:satz:geometrischeinfach} garantiert, +dass der Eigenwert einfach ist. +Es ist aber immer noch möglich, dass die algebraische Vielfachheit +von $\varrho(A) >1$ ist, dass also $\dim\mathcal{E}_{\varrho(A)}(A)>1$ +ist. +Dies ist jedoch nicht der Fall. + \begin{satz} \label{buch:wahrscheinlichkeit:satz:algebraischeinfach} -Der verallgemeinerte Eigenraum zum Eigenwert $\varrho(A)$ einer -positiven Matrix $A$ ist eindimensional. +Sei $A$ eine positive Matrix und $p^t$ ein positiver Eigenvektor +der Matrix $A^t$ zum Eigenwert $\varrho(A^t)=\varrho(A)$. Ist $u$ der Eigenvektor von $A$ zum Eigenwert $\varrho(A)$ nach -Satz~\ref{buch:wahrscheinlichkeit:satz:geometrischeinfach} -und $p^t$ der entsprechende Eigenvektor $A^t$, dann -ist +Satz~\ref{buch:wahrscheinlichkeit:satz:geometrischeinfach}, +dann ist \[ \mathbb{R}^n = @@ -643,6 +691,8 @@ ist \ker p \] eine Zerlegung in invariante Unterräume von $A$. +Insbesondere ist der verallgemeinerte Eigenraum $\mathcal{E}_{\varrho(A)}(A)$ +von $A$ eindimensional. \end{satz} \begin{proof}[Beweis] @@ -652,7 +702,8 @@ Insbesondere ist $u\not\in\ker p$ Es ist klar, dass $A\langle u\rangle = \langle Au\rangle = \langle u\rangle$ ein invarianter Unterraum ist. -Für einen Vektor $x\in\mathbb{R}^n$ mit $px=0$ erfüllt das Bild $Ax$ +Für einen Vektor $x\in\mathbb{R}^n$ mit $px=0$, also $x\in\ker p$, +erfüllt das Bild $Ax$ die Gleichung \[ p(Ax)=(pA)x=(A^tp^t)^tx= \varrho(A)(p^t)^tx @@ -666,8 +717,8 @@ $\ker p$ ist $(n-1)$-dimensional, $\langle u\rangle$ ist eindimensional und $u$ ist nicht in $\ker p$ enthalten. Folglich spannen $\langle u\rangle$ und $\ker p$ den ganzen Raum auf. -Gäbe es einen weitern linear unabhängigen Vektor im verallgemeinerten -Eigenraum von $\mathcal{E}_{\varrho(A)}$, dann müsste es auch einen +Gäbe es einen weiteren linear unabhängigen Vektor im verallgemeinerten +Eigenraum $\mathcal{E}_{\varrho(A)}(A)$, dann müsste es auch einen solchen Vektor in $\ker p$ geben. Da $\ker p$ invariant ist, müsste es also auch einen weiteren Eigenvektor $u_2$ zum Eigenwert $\varrho(A)$ in $\ker p$ geben. @@ -712,10 +763,10 @@ Dann ist $\varrho(A)$ der einzige Eigenwert vom Betrag $\varrho(A)$ und er hat geometrische und algebraische Vielfachheit $1$. \end{satz} -\begin{proof}[Beweis] +\begin{proof}[Beweisansatz] Nach Voraussetzung gibt es ein $n$ derart, dass $A^n>0$. Für $A^n$ gelten die Resultate von Satz~\ref{buch:wahrscheinlichkeit:satz:perron-frobenius}. - -XXX TODO +Man kann zeigen, dass die Eigenvektoren von $A^n$ auch +Eigenvektoren von $A$ sind. \end{proof} diff --git a/buch/chapters/90-crypto/chapter.tex b/buch/chapters/90-crypto/chapter.tex index d2fcbbf..2ea0932 100644 --- a/buch/chapters/90-crypto/chapter.tex +++ b/buch/chapters/90-crypto/chapter.tex @@ -5,14 +5,14 @@ % (c) 2020 Prof Dr Andreas Müller, Hochschule Rapperswil % % !TeX spellcheck = de_CH -\chapter{Anwendungen in Kryptographie und Codierungstheorie +\chapter{Anwendungen in Kryptographie \label{buch:chapter:kryptographie}} -\lhead{Kryptographie und Codierungstheorie} +\lhead{Kryptographie} \rhead{} -Die algebraische Theorie der endlichen Körper hat sich als besonders -nützliche herausgestellt in der Krypographie. +Die algebraische Theorie der endlichen Körper hat sich +in der Krypographie als besonders nützliche herausgestellt. Die Eigenschaften dieser Körper sind reichhaltig genug, um -kryptographsch widerstandsfähige Algorithmen zu liefern, die +kryptographisch widerstandsfähige Algorithmen zu liefern, die auch in ihrer Stärke beliebig skaliert werden können. Gleichzeitig liefert die Algebra auch eine effiziente Implementierung. In diesem Abschnitt soll dies an einigen Beispielen gezeigt werden. diff --git a/buch/chapters/references.bib b/buch/chapters/references.bib index fb88d09..979f985 100644 --- a/buch/chapters/references.bib +++ b/buch/chapters/references.bib @@ -162,3 +162,35 @@ abstract = "In this paper, we present Google, a prototype of a large-scale searc year = 2011, language = {english} } + +@article{francis:watkins_paper, + ISSN = {00029890, 19300972}, + URL = {http://www.jstor.org/stable/10.4169/amer.math.monthly.118.05.387}, + author = {David S. Watkins}, + journal = {The American Mathematical Monthly}, + number = {5}, + pages = {387--403}, + publisher = {Mathematical Association of America}, + title = {Francis’ Algorithm}, + volume = {118}, + year = {2011} +} + +@article{skript:pftopo, + author = { Geon H. Ghoe }, + journal = { Communications of the Korean Mathematical Society }, + volume = {9}, + number = {3}, + pages = {565--570}, + year = {1994} +} + +@article{skript:hilbertmetric, + author = { Elon Kohlberg and John W. Pratt }, + title = {The contraction mapping approach to the Perron-Frobenius theory: why Hilbert's metric? }, + journal = {Mathematics of Operations Research}, + year = {1982}, + volume = {7}, + number = {2}, + pages = {192--210} +} diff --git a/buch/papers/erdbeben/Teil_Fabio.tex b/buch/papers/erdbeben/Teil_Fabio.tex index 653470a..b51d1f0 100644 --- a/buch/papers/erdbeben/Teil_Fabio.tex +++ b/buch/papers/erdbeben/Teil_Fabio.tex @@ -1,106 +1,79 @@ \section{Anwendung des Kalman-Filters} -\subsection{Ziel} -Bis jetzt haben wir gelesen, was das Kalman-Filter bewirkt und wie es funktioniert. -Nun möchten wir mit einem Beispiel herausfinden, ob das Filter unsere gesuchte Grösse $f(t)$ bestimmen kann. +Bis jetzt haben wir gesehen, was das Kalman-Filter bewirkt und wie es funktioniert. +Nun möchten wir mit einem konkreten Beispiel herausfinden, +ob das Filter unsere gesuchte Grösse $f(t)$ bestimmen kann. -\subsection{Künstliche Erdbebendaten} -Da wir keine Rohdaten über vergangene Erdbeben zur Hand haben, müssen wir mittels Matlab künstliche Daten erzeugen und sie dann in das Filter eingeben. +Da wir keine Rohdaten über vergangene Erdbeben zur Hand haben, +müssen wir mittels Simulation künstliche Daten erzeugen. +Diese können wir dann mit unserem Filter verarbeiten. Diese Vorgehensweise erlaubt uns das Erdbeben beliebig zu gestalten -und weil es digital simuliert wird, haben wir keine Bauschäden zu beklagen. +und weil es digital simuliert wird, haben wir auch keine Bauschäden zu beklagen. \subsection{Wahl der Schwingung} Wir müssen uns überlegen, mit welcher Schwingung wir ein realitätsnahes Beben erzeugen können. - Mit einer ungedämpften harmonischen Schwingung können wir zwar die meisten Vorgänge in der Physik erklären. -Da aber unser Erdbeben irgendwann abklingen muss, wählen wir die gedämpfte harmonische Schwingung. -Die dazugehörige Schwingungsgleichung lautet - -\begin{equation} - y = A e^{-\lambda t} \sin(\omega t) -\end{equation} - -Für die Variablen der harmonisch gedämpften Schwingung setzen wir die Werte - -\begin{equation} -A = 5 -\end{equation} - -ein. - -$A$ ist die Amplitude der Schwingung, die uns die Heftigkeit des Erdebebens beschreibt. -Sie ist vergleichbar mit der Magnitude. - -$\omega$ definiert sich durch - +Da ein Erdbeben vorteilhafterweise irgendwann abklingen sollte, +wählen wir eine gedämpfte harmonische Schwingung \begin{equation} - \omega = 2 \pi f + y = A e^{-\lambda t} \sin(2\pi f t). \end{equation} -wobei die Frequenz $f$ mit +In unsere Simulation können wir die Parameter frei wählen. +Wir setzten $A = 5$ für die Amplitude der Schwingung. +Sie beschreibt die Heftigkeit des Erdebebens und ist vergleichbar mit der Magnitude. +Für die Frequenz $f$ wählen wir eine Zufalls-Sequenz mit Erwartungswert und Standardabweichung \begin{equation} - f = E(\mathrm{Frequenz}) + \sigma^2(\mathrm{Frequenz}) + \mu = \SI{15}{\hertz} + \qquad\text{und}\qquad + \sigma = \SI{10}{\hertz}. \end{equation} - -erzeugt wird. - -Zusätzlich haben wir $f$ mit dem Savitzky-Golay-Filter gefiltert. -Das Savitzky-Golay-Filter schaut sich immer eine definierte Anzahl von Datenpunkte an +Zusätzlich haben wir $f$ mit einem Savitzky-Golay-Filter gefiltert. +Ein Savitzky-Golay-Filter schaut sich immer eine definierte Anzahl von Datenpunkte an und bildet ein Polynom $n$-ter Ordnung. -In unserer Anwendung schaut sich das Filter, im Sinne eines verschieblichen Fensters, -jeweils zehn aufeinanderfolgende Datenpunkte an und bildet ein Polynom $0$-ter Ordnung. -Da wir den Grad $0$ gewählt haben, erhalten wir pro zehn Punkte eine Gerade mit der Steigung $0$. -Diese Art von der Filterung nennt sich gleitender Mittelwert. - -Für den Erwartungswert und die Standardabweichung setzen wir die Zahlen - -\begin{equation} -E(f) = \SI{15}{\hertz} -\end{equation} - -und -\begin{equation} -\sigma^2 = \SI{10}{\hertz} -\end{equation} - -ein. +In unserer Anwendung schaut sich das Filter, im Sinne eines verschiebbaren Fensters, +jeweils elf aufeinanderfolgende Datenpunkte an und bildet ein Polynom $0$-ter Ordnung, +also eine Konstante. +Somit erhalten wir mit Matlab-Standardfunktionen einen gleitenden Mittelwert, +um all zu schnelle Änderungen der Frequenz zu unterdrücken. $\lambda$ ist die Bodendämpfung, für die wir $0.2$ wählen. -Sie ist dafür verantwortlich, dass unser Erdbeben abklingen wird und kreiert bei der gedämpften Schwingung die typische Hüllkurve der Amplitude. +Sie ist dafür verantwortlich, dass unser Erdbeben abklingt +und kreiert bei der gedämpften Schwingung die typische Hüllkurve. Wir nehmen an, dass $\lambda$ ein Materialparameter von geologischen Böden ist. \subsection{Versuch im Standardfall} Im nächsten Schritt müssen wir sinnvolle Systemparameter für unseren Seismographen definieren. Eine kurze Recherche zeigt, dass die Masse ein Gewicht von ca.\ \SI{100}{\gram} hat. -Zur Federkonstante D und Dämpfung k konnten wir leider keine brauchbaren Grössen finden und treffen die Annahme, dass $D = 1$ und $k = 0.01$. +Zur Federkonstante $D$ und Dämpfung $k$ konnten wir leider keine brauchbaren Grössen finden. +Wir treffen die Annahmen $D = 1$ und $k = 0.01$. Für die Masse definieren wir $m = 0.01$. Für das Prozessrauschen werden die Bedingungen - \begin{equation} - Q = \left( - \begin{array}{ccc} - {\sigma_x }^2& 0& 0 \\ - 0 & {\sigma_v }^2& 0\\ - 0 & 0& {\sigma_f }^2\\ - \end{array}\right)= \left( - \begin{array}{ccc} - {0.00001}^2& 0& 0 \\ - 0 & {0.00001}^2& 0\\ - 0 & 0& {1 }^2\\ - \end{array}\right) + Q = + \begin{pmatrix} + \sigma_x ^2 & 0 & 0 \\ + 0 & \sigma_v^2 & 0\\ + 0 & 0 & \sigma_f^2 \\ + \end{pmatrix}= + \begin{pmatrix} + 0.00001^2& 0& 0 \\ + 0 & 0.00001^2& 0\\ + 0 & 0& 1^2 \\ + \end{pmatrix} \end{equation} - angesetzt. - +Die Annahme, dass sich die Erdbebenkraft $f$ nicht ändert, +kompensieren wir hier endlich durch einen grossen Wert von $\sigma_f^2$. Auch für die Messung setzen wir ein Rauschen voraus und definieren - \begin{equation} -R= ({\sigma_x}^2)= -({0.00001}^2) +R= (\sigma_x^2)= +(0.00001^2). \end{equation} - -Sind nun die benötigten Systemparameter und das Rauschen definiert, erzeugen wir das Erdbeben und schauen, wie gut das Kalman-Filter die äussere Beschleunigung schätzen kann. +Damit sind nun die benötigten Systemparameter und das Rauschen definiert. +Als nächstes erzeugen wir ein Erdbeben und schauen, +wie gut das Kalman-Filter die äussere Beschleunigung schätzen kann. \subsection*{Ergebnis} @@ -112,18 +85,25 @@ Zoomen wir näher ran, erkennen wir wieder im Positions-Diagramm eine Überlager Die Masse schwingt mit einer tiefer Frequenz und hoher Amplitude, hingegen das Erdbeben mit einer hohen Frequenz und tiefer Amplitude. Vergleichen wir nun die Position mit der Kraft, stellen wir fest, dass das Kalman-Filter eine Schätzung wiedergibt, die auch eine Frequenz von \SI{15}{\hertz} hat. -Das Filter war imstande die Eigenfrequenz zu eliminieren und die tatsächliche Kraft des Erdbebens zu wiedergeben. +Das Filter war imstande die Eigenfrequenz zu eliminieren und die tatsächliche Kraft des Erdbebens wiederzugeben. \begin{figure} - \begin{center} - \includegraphics[width=\linewidth,keepaspectratio]{papers/erdbeben/images/standard.PDF} - \caption{Das Position-Zeit-Diagramm zeigt uns die typische Aufzeichnung eines Seismographen während eines Erdbebens. Um die Geschwindigkeit zu erhalten müssen wir die Position einmal ableiten. Ein weiteres Ableiten erzeugt uns die Beschleunigung, respektive die Kraft. Sehr gut ersichtlich ist die Hüllkurve der Amplitude, wie wir sie bei einer gedämpften Schwingung erwarten. Erst das Vergrössern an die Datenpunkte zeigt uns auf, wie gut die Schätzung des Kalman-Filters funktioniert.} + \begin{center} + \includegraphics[width=.95\linewidth,keepaspectratio]{papers/erdbeben/images/standard.PDF} + \caption{ + Das Position-Zeit-Diagramm zeigt eine typische Aufzeichnung eines Seismographen während eines Erdbebens. + Sehr gut ersichtlich ist die Hüllkurve, wie wir sie bei einer gedämpften Schwingung erwarten. + In der Vergrösserung wird die Überlagerung aus Eigenschwingung und Erdbeben gut ersichtlich. + Die Geschwindigkeit und schliesslich die Kraft weden aus der Position durch unser Kalman-Filter geschätzt. + Erst das Vergrössern an die Datenpunkte zeigt, wie gut die Schätzung des Kalman-Filters funktioniert. + In der Kraft ist die Eigendynamim nicht mehr ersichtlich. Unser Filter funktioniert. + } \label{erdbeben:fig:standard-alles} - \end{center} + \end{center} \end{figure} \subsection{Veränderung der Systemparameter} -Was wir nun austesten möchten, sind die Auswirkungen wenn z.B. der Seismograph andere Systemparameter aufweist. +Wir möchten nun testen, was die Auswirkungen sind, wenn zum Beispiel der Seismograph andere Systemparameter aufweist. Wir nehmen an, dass sich im Vergleich zum Standardfall die Masse erhöht, die Federkonstante schwächer und die Bodendämpfung doppelt so stark wirkt. Somit gilt neu \[ @@ -133,67 +113,97 @@ D = 0.5 \qquad \text{und} \qquad k = 0.02. \] +Da wir mit dieser Anpassung die Trägheit des Seismogrammes erhöht haben, +erwarten wir eine langsamere Bewegung der Masse, das heisst die Eigenfrequenz wird reduziert. -Da wir mit dieser Anpassung die Trägheit des Seismogrammes erhöht haben, erwarten wir sicher eine langsamere Bewegung der Masse, das heisst die Frequenz wird sich reduzieren. - -Betrachten wir die Abbildung~\ref{erdbeben:fig:systemparameter-geaendert} können wir diese Erwartung bestätigen. +Betrachten wir Abbildung~\ref{erdbeben:fig:systemparameter-geaendert} können wir diese Erwartung bestätigen. Nebst dem bemerken wir eine grössere Auslenkung der Position, die wir auf die höhere Energie der Masse und geringeren Rücklenkkraft der Feder begründen können. \begin{figure} - \begin{center} - \includegraphics[width=\linewidth,keepaspectratio]{papers/erdbeben/Systemparameter_geaendert_2.PNG} - \caption{Im Geschwindigkeits-Diagramm erkennen wir, dass sich im Vergleich zum Standardfall, die Auslenkung und Frequenz vergrössert hat. Dies wird mit der Erhöhung der Masse und somit der Trägheit begründet. Auch stellen wir fest, dass die Positionsmessung überwiegend die Eigenfrequenz misst.} + \begin{center} + \includegraphics[width=\linewidth,keepaspectratio]{papers/erdbeben/Systemparameter_geaendert.PNG} + \caption{ + Im Geschwindigkeits-Diagramm erkennen wir, + dass sich im Vergleich zum Standardfall die Auslenkung und Frequenz vergrössert hat. + Dies wird mit der Erhöhung der Masse und somit der Trägheit begründet. + Auch stellen wir fest, dass die Positionsmessung überwiegend die Eigenfrequenz misst. + Diese wird in der Schätzung der Kraft dennoch sehr gut kompensiert. + } \label{erdbeben:fig:systemparameter-geaendert} - \end{center} + \end{center} \end{figure} -<<<<<<< Updated upstream -======= - - - ->>>>>>> Stashed changes \subsection{Verstärkung des Prozessrauschens} Falls wir unseren Seismographen in der Nähe einer grösseren Stadt aufstellen, so müssen wir aufgrund der Vibrationen mit einem stärkeren Prozessrauschen rechnen. Dieses Rauschen beeinflusst die Varianzen der Position und Geschwindigkeit in der Matrix $Q$. -\ref{erdbeben:Wahrscheinlichkeit} -Aus diesem Grund erhöhen wir die Standardabweichungen in der Matrix $Q$ um den Faktor $100$. -Die Auswertung in Abbildung~\ref{erdbeben:fig:prozessrauschen-geaendert} zeigt auf, dass das Kalman-Filter die Schätzung der Kraft nur gering an den Messwerten anpasst. -Die Theorie dazu haben wir im Kapitel \nameref{Wahrscheinlichkeit} angeschaut. +Aus diesem Grund erhöhen wir die Standardabweichungen der Positions $\sigma_s$ und Geschwindigkeit $\sigma_v$ in der Matrix $Q$ um den Faktor $100$. +Die Auswertung in Abbildung~\ref{erdbeben:fig:prozessrauschen-geaendert} zeigt auf, +dass das Kalman-Filter die Schätzung der Kraft nur gering an den Messwerten anpasst, +da wir den Schätzungen für die Position nun wenig Vertrauen schenken und stärker der Modell-Annahme $\dot f = 0$ folgen. +Die Theorie dazu haben wir im Abschnitt~\ref{erdbeben:Wahrscheindlichkeit} angeschaut. \begin{figure} - \begin{center} - \includegraphics[width=\linewidth,keepaspectratio]{papers/erdbeben/images/Prozessrauschen_geaendert.PDF} - \caption{Mit dem Erhöhen des Prozessrauschens gehen wir von einer grösseren Unsicherheit der Systemmatrix aus. Aus diesem Grund folgt das Filter vor allem den Messwerten, was sichtbare Folgen für die Schätzkurve im Kraft-Zeit-Diagramm hat. Hier möchte das Filter auch den Messwerten folgen. Da wir aber für die Kraft keine Messwerte aufzeichnen, erhalten wir eine sehr schwache Kurve. Die Position kann immernoch präzise geschätzt werden und die Ableitung zur Geschwindigkeit ergibt gute Resultate. Jedoch ist die Schätzkurve der Kraft sehr weit von der idealen Kurve entfernt und nicht nutzbar.} + \begin{center} + \includegraphics[width=.95\linewidth,keepaspectratio]{papers/erdbeben/images/Prozessrauschen_geaendert.PDF} + \caption{ + Mit dem Erhöhen des Prozessrauschens gehen wir von einer grösseren Unsicherheit der Systemmatrix aus. + Aus diesem Grund folgt das Filter vor allem den Messwerten, + was sichtbare Folgen für die Schätzkurve im Kraft-Zeit-Diagramm hat. + Hier möchte das Filter auch den Messwerten folgen. + Da wir aber für die Kraft keine Messwerte aufzeichnen, + erhalten wir eine sehr schwache Kurve. + Die Position kann immernoch präzise geschätzt werden und die Ableitung zur Geschwindigkeit ergibt gute Resultate. + Jedoch ist die Schätzkurve der Kraft sehr weit von der idealen Kurve entfernt und nicht nutzbar. + } \label{erdbeben:fig:prozessrauschen-geaendert} - \end{center} + \end{center} \end{figure} \subsection{Verstärkung des Messrauschens} -\ref{erdbeben:Wahrscheinlichkeit} Als letztes verstärken wir das Messrauschen um den Faktor $100$ und belassen wieder den Rest wie im Standardfall. -Wie man eigentlich schon erwarten kann, zeigt uns die Abbildung~\ref{erdbeben:fig:messrauschen-geaendert}, dass das Signal des Messsensors vom Messrauschen gestört wird. -Weil die Messung somit ungenau wird, kann das Kalman-Filter nicht mehr genau arbeiten und produziert einen ungenauen Output. -Auch hier verweisen wir auf Kapitel \nameref{Wahrscheinlichkeit}. +Wie man eigentlich schon erwarten kann, +zeigt uns die Abbildung~\ref{erdbeben:fig:messrauschen-geaendert}, +dass das Signal des Messsensors vom Messrauschen start gestört wird. +Weil die Messung zu ungenau ist, +kann das Kalman-Filter nicht mehr gut arbeiten und produziert einen ungenauen Output. \begin{figure} - \begin{center} - \includegraphics[width=\linewidth,keepaspectratio]{papers/erdbeben/images/Messrauschen_geaendert.PDF} - \caption{Im Kraft-Zeit-Diagramm erhalten wir nur bis ca. $t = 10$ gute Schätzwerte. Von $t = 10$ bis $t = 30$ wirkt das Messrauschen zu stark und erhalten keine brauchbaren Werte mehr. Im Position-Zeit-Diagramm erhielten wir bis jetzt immer genaue Schätzungen. Mit einem starken Messrauschen fällt es nun dem Filter schwerer, präzise Werte zu generieren. Die Nahaufnahme im Kraft-Zeit-Diagramm bestätigt uns aber, dass die Messfehler zu gross sind, um ein klares Bild über die äussere Kraft zu erhalten.} - \label{erdbeben:fig:messrauschen-geaendert} - \end{center} + \begin{center} + \includegraphics[width=.95\linewidth,keepaspectratio]{papers/erdbeben/images/Messrauschen_geaendert.PDF} + \caption{ + Im Kraft-Zeit-Diagramm erhalten wir nur bis ca.\ $t = 10$ gute Schätzwerte. + Ab $t = 10$ wirkt das Messrauschen zu stark und wir erhalten keine brauchbaren Werte mehr. + Im Position-Zeit-Diagramm erhielten wir bis jetzt immer genaue Schätzungen. + Mit einem starken Messrauschen fällt es dem Filter nun schwerer, + präzise Schätzungen zu berechnen. + Die Nahaufnahme im Kraft-Zeit-Diagramm bestätigt uns, + dass die Messfehler zu gross sind, + um ein klares Bild über die äussere Kraft zu erhalten. + } + \label{erdbeben:fig:messrauschen-geaendert} + \end{center} \end{figure} \subsection{Zusammenfassung} -Wir haben uns zum Ziel gesetzt, die äussere Beschleunigung $a(t)$, bzw. die Kraft $f(t)$ eines Erdbebens zu ermitteln. +Wir haben uns zum Ziel gesetzt, +die äussere Beschleunigung $a(t)$, +beziehungsweise die Kraft $f(t)$ eines Erdbebens +aus den Messugnen eines Seismographen zu berechen. + +Wir haben einen Seismographen mathematisch beschrieben und +mit der Software Matlab Messresultate während eines künstlichen Erdbebens erzeugt. +Diese Messwerte haben wir mit einem Kalman-Filter bearbeitet, +um aus den Messwerten wieder das Erdbeben zu gewinnen. -Mit der Software Matlab haben wir einen virtuellen Seismographen gebaut und ein künstliches Erdbeben erzeugt. -Der Seismograph war fähig die Position der Masse während der Einwirkung des Erdbebens aufzuzeichnen. +Der Seismograph war fähig, die Position der Masse während der Einwirkung des Erdbebens aufzuzeichnen. $a(t)$ kann zwar nicht mit Sensoren gemessen werden, jedoch erhalten wir $a(t)$ durch zweifaches Ableiten. Da wir so aber die innere Beschleunigung erhalten, mussten wir das Kalman-Filter anwenden. -Das Kalman-Filter half uns die äussere Beschleunigung zu schätzen und lieferte erstaunlich genaue Werte. +Das Kalman-Filter half uns, die äussere Beschleunigung zu schätzen, und lieferte erstaunlich genaue Werte. Ausserdem hat es das Filter geschafft, die Eigenfrequenz der Masse und die Erdbebenfrequenz zu separieren. Folglich erhielten wir eine Schätzung, die nur das Erdbeben betraf. -Zuletzt haben wir aufgezeigt, das Veränderungen an den System- und Rauschparametern die Genauigkeit und Zuverlässigkeit des Kalman-Filters beeinträchtigen können. +Zuletzt haben wir aufgezeigt, +dass Veränderungen an den System- und Rauschparametern die Genauigkeit und Zuverlässigkeit +des Kalman-Filters beeinträchtigen können. +Wir haben gesehen, dass aus zu schlechten Sensordaten auch mittels Filterung keine genauen Aussagen möglich sind. diff --git a/buch/papers/erdbeben/references.bib b/buch/papers/erdbeben/references.bib index 444c82d..9bcf57d 100644 --- a/buch/papers/erdbeben/references.bib +++ b/buch/papers/erdbeben/references.bib @@ -28,6 +28,13 @@ title = {Kalmanfilter}, urldate = {2021-07-0}} +@article{erdbeben:skript:wrstat, + author = {Andreas Müller}, + url = {https://github.com/AndreasFMueller/WrStat}, + date = {2021-09-11}, + title = {Wahrscheinlichkeitsrechnung und Statistik} +} + @online{erdbeben:bibtex, date = {2020-02-06}, day = {6}, diff --git a/buch/papers/erdbeben/teil0.tex b/buch/papers/erdbeben/teil0.tex index d32b316..a3fa6a5 100644 --- a/buch/papers/erdbeben/teil0.tex +++ b/buch/papers/erdbeben/teil0.tex @@ -3,11 +3,9 @@ % % (c) 2020 Prof Dr Andreas Müller, Hochschule Rapperswil %% -\section{Was ist ein Erdbeben? \label{erdbeben:section:teil0}} -\rhead{Erdbeben} -Für das Verständnis möchten wir zuerst erklären, was ein Erdbeben genau ist. -Das soll uns helfen, eine Verknüpfung zwischen dem Naturphänomen und der mathematischen Problemstellung herzustellen. +\rhead{Erdbeben} +\noindent Unter einem Erdbeben verstehen wir eine Erschütterung des Erdkörpers. Dabei reiben zwei tektonische Platten aneinander, welche sich durch die Gesteinsverzahnung gegenseitig blockieren. Diese Haftreibung durch die Steine wird so lange aufgebaut, bis sie nicht mehr gehalten werden kann. @@ -15,17 +13,31 @@ Wenn dies passiert, entlädt sich die aufgebaute Spannung und setzt enorme Energ Ein Erdbeben breitet sich vom Erdbebenherd in allen Richtungen gleich aus. Vergleichbar ist, wenn man einen Stein in einen Teich wirft und die Wellen beobachten kann, die sich ausbreiten. -\subsection{Funktion eines Seismograph} +\section{Funktion eines Seismographen} Um ein Erdbeben kenntlich zu machen, werden in der Regel Seismographen mit vielen Sensoren verwendet. -Ein Seismograph besteht im Grunde aus einer federgelagerten Masse. Wirkt eine Bodenerregung auf das Gerät ein, schwing das Gehäuse und dadurch auch die gekoppelte Masse. -Stoppt das Erdbeben, schwingt das Gehäuse nicht mehr. -Die Masse schwing jedoch in seiner Eigendynamik weiter. -Eine Relativbewegung des Bodens kann damit als Auslenkung im Zeitverlauf gemessen werden. -In modernen Seismographen wird die Bodenbewegung in alle Richtungen gemessen, sowohl Horizontal als auch Vertikal. -Wir konstruieren uns eine einfachere Version eines Seismographen mit einem Gehäuse, an dem zwei Federn und eine Masse befestigt sind. -Der Seismograph ist in Abbildung ~\ref{erdbeben:Seismograph} ersichtlich. -Ein Sensor unter der Masse misst die Position, bzw. die Auslenkung der Feder und der Masse. -Dies bedeutet, unser Seismograph kann nur in eine Dimension Messwerte aufnehmen. +Ein Seismograph besteht im Grunde aus einer federgelagerten Masse. +Bei einem Erdbeben folgt das Gehäuse direkt der Bewegung des Erdbebens. +Die federgelagerte Masse wird jedoch erst durch die Feder bewegt und folgt verzögert. +Zudem schwingt die Masse auch ohne Erdbeben weiter -- das System besitzt eine Eigendynamik. +Eine Relativbewegung des Gehäuses kann folglich als Auslenkung im Zeitverlauf gemessen werden. +Allerdings misst man so nicht direkt das Erbeben, sondern eine Überlagerung der Effekte aus Erdbeben- und Federkraft. + +In modernen Seismographen wird die Bodenbewegung in alle Richtungen gemessen, +sowohl horizontal als auch vertikal. +Wir konstruieren hier eine einfachere Version eines Seismographen mit einem Gehäuse, +an dem zwei Federn und eine Masse befestigt sind. +Abbildung~\ref{erdbeben:Seismograph} zeigt eine schematische Darstellung unseres Systems. +Ein Sensor unter der Masse misst die Position der Masse relativ zum Gehäuse. +Unser Seismograph misst also nur eindimensional. + +Für mehrere Dimensionen würde der Satz von Pythagoras für die Auslenkung der Federn benötigt. +Die benötigten Quadrate und Wurzeln brechen jedoch die Linearität des Systems. +Die Systembeschreibung wird dann deutlich komplexer, bringt aber nichts wesentlich Neues hervor. +Wir beschränken uns deshalb auf den linearen Fall. + +Wir werden sehen, dass diese Art der Problemstellung effektiv mittels Kalman-Filter gelöst werden kann. +Für ein nicht-lineares System werden Extended Kalman-Filter benötigt, +bei denen die System-Matrix $A$ durch die Jacobi-Matrix ersetzt wird. \begin{figure} \begin{center} @@ -35,75 +47,69 @@ Dies bedeutet, unser Seismograph kann nur in eine Dimension Messwerte aufnehmen. \end{center} \end{figure} -\subsection{Ziel} -Unser Seismograph misst nur die Position der Masse über die Zeit. -Wir wollen jedoch die Beschleunigung $a(t)$ des Boden, bzw. die Kraft $f(t)$, welche auf das Gehäuse wirkt, bestimmten. -Anhand dieser Beschleunigung, bzw. der Krafteinwirkung durch die Bodenbewegung, wird später das Bauwerk bemessen. +Unser Seismograph misst jedoch nur die Position der Masse über die Zeit. +Wir wollen aber die Beschleunigung $a(t)$ des Boden, +respektive die Kraft, +welche auf das Gehäuse wirkt, bestimmen. +Anhand dieser Beschleunigung, +beziehungsweise der Krafteinwirkung durch die Bodenbewegung, +wird später das Bauwerk bemessen. Dies bedeutet, die für uns interessante Grösse $f(t)$ wird nicht durch einen Sensor erfasst. Jedoch können wir durch zweifaches ableiten der Positionsmessung $s(t)$ die Beschleunigung der Masse berechnen. -Das heisst: Die Messung ist zweifach Integriert die Kraft $f(t)$ inklusive der Eigendynamik der Masse. -Um die Krafteinwirkung der Masse zu berechnen, müssen wir Gleichungen für unser System finden. +Die Messung entspricht also dem zweiten Integral der Kraft $f(t)$, +wobei diese einerseits durch das Erdbeben, und andererseits durch die Federn zustande kommt. +Im Folgenden möchten wir die Erdbeben- und Federkräfte trennen. +Dafür benötigen wir zuerst eine mathematische Beschreibung unseres Systems. \subsection{Systemgleichung} Im Paper~\cite{erdbeben:mendezmueller} wurde das System gleich definiert und vorgegangen. Im Fall unseres Seismographen, handelt es sich um ein Feder-Masse-Pendel. -Dieser kann durch die Differentialgleichung zweiter Ordnung einer gedämpften Schwingung am harmonischen Oszillator beschrieben werden. -Die Gleichung lautet: +Dieses kann als gedämpfter harmonischer Oszillator beschrieben werden. +Die zugehörige Differentialgleichung lautet: \begin{equation} + \label{erdbeben:Systemgleichung} m\ddot s + 2k \dot s + Ds = f. \end{equation} wobei $m$ die Masse, $k$ die Dämpfungskonstante und $D$ die Federkonstante bezeichnet. - -Da die Differentialgleichung linear ist möchten wir diese Gleichung in die Darstellung $\dot x = Ax$ überführen, wobei $x$ der Zustandsvektor und $A$ die Systemmatrix bezeichnet. Dazu verwenden wir die Subsitution: -\[ -s_1 = s -\qquad \text{und} \qquad -s_2 = \dot s. -\] -Somit entstehen die Gleichungen für die Geschwindigkeit $ \dot s_1(t)$ der Masse : -\[ \dot {s_1} = {s_2}\] -und -\[ \dot s_2 = -\frac{D}{m} {s_1} -\frac{2k}{m} {s_2} + \frac{f} {m} \] -für die Beschleunigung $\dot s_2(t)$ der Masse. -Diese können wir nun in der Form -\[ \ddot f =-\frac{D}{m} {s_1} -\frac{2k}{m} {s_2} + \frac{f} {m} \] -als skalare Gleichung darstellen. - -Die für uns relevanten Zustände sind die Position der Masse, die Geschwindigkeit der Masse und die äussere Beschleunigung des ganzen Systems. -Unüblich ist nun, dass der Stör-Term $f$ in Gleichung (20.1) gerade das ist, was wir eigentlich bestimmen möchten. -In unserem Fall wird die äusseren Beschleunigung gesucht, da diese der Erdbebenanregung gleich kommt. +Für lineare Systeme ist eine Matrix-Darstellung handlicher. +Wir möchten diese Gleichung folglich in die Darstellung $\dot x = Ax$ überführen, +wobei $x$ der Zustandsvektor und $A$ die Systemmatrix bezeichnet. +Wir subsituieren $\dot s = v$ für die Geschwindigkeit und erhalten das Gleichungssystem +\begin{align} + \begin{split} + \dot s &= v \\ + \dot v &= -\frac{D}{m} {s} -\frac{2k}{m} {v} + \frac{f} {m}. + \end{split} + \label{erdbenen:systemgleichungen} +\end{align} + +Die relevanten Zustände sind also die Position $s$ und die Geschwindigkeit $v$. +Die für uns eigentlich interessante Grösse ist jedoch der Stör-Term $f$. +Dieser entspricht der Kraft durch das Erdbeben. Deshalb nehmen wir $f$ als dritte Grösse in den Zustandsvektor auf und definieren: - \[ - x = (s_1, s_2, f)^T. + x = \begin{pmatrix} {s} \\ {v} \\ {f} \end{pmatrix} \] -Für die Standard-Form $\dot x = Ax$ brauchen wir als nächstes die Ableitungen aller Elemente von $x$. Für $\dot s_1$ und $\dot s_2$ folgen diese direkt aus Gleichung (20.1), aber über $\dot f$ wissen wir nichts. -Wir müssen also eine Annahme treffen: $\dot f = 0$. Diese Annahme ist im Allgemeinen falsch, aber etwas Besseres haben wir zurzeit nicht zur Verfügung. -Zudem treffen wir die Annahme, das sich die Kraft über die Beobachtungszeit nicht verändert. +Für die Standard-Form $\dot x = Ax$ brauchen wir als nächstes die Ableitungen aller Elemente von $x$. +Für $s$ und $v$ haben wir diese in Gleichung~\eqref{erdbenen:systemgleichungen} bereits gefunden. +Über die Kraft $f$ wissen wir jedoch nichts. +Wir müssen also eine Annahme treffen: Die Kraft ändert sich nicht, $\dot f = 0$. +Diese Annahme ist im Allgemeinen natürlich falsch, aber etwas Besseres haben wir nicht zur Verfügung. Wir werden dies in einem späteren Schritt kompensieren müssen. -Da die Kraft unbekannt ist, wird die letzte Zeile mit Nullen gefüllt, denn genau diese Werte wollen wir. - -Durch Rücksubstituion ergibt sich uns folgende Systemgleichung in Matrix schreibweise, wobei $\dot {s_1}= v$ ist. Damit haben wir nun alles, was wir für die Matrix-Darstellung von Gleichung (20.1) benötigen. Diese lautet: +Wir haben nun alles für die Matrix-Form von Gleichung~\eqref{erdbeben:Systemgleichung} zusammen. +Sie lautet: \begin{equation} -\frac{d}{dt} \left(\begin{array}{c} s(t) \\ v(t) \\ f(t) \end{array}\right) = \left( - \begin{array}{ccc} -0 & 1& 0 \\ -- \frac{D}{m} &-\frac{2k}{m} & \frac{1} {m}\\ -0 & 0 & 0\\ -\end{array}\right) \left(\begin{array}{c} s(t)\\ v(t)\\ f(t) \end{array}\right). + \frac{d}{dt} \begin{pmatrix} s(t) \\ v(t) \\ f(t) \end{pmatrix} + = + \begin{pmatrix} + \phantom- 0 & \phantom-1& 0 \\ + - \frac{D}{m} &-\frac{2k}{m} & \frac{1} {m} \\ + \phantom-0 & \phantom-0 & 0\\ + \end{pmatrix} + \begin{pmatrix} s(t) \\ v(t) \\ f(t) \end{pmatrix}. + \label{erdbeben:systemmatrix} \end{equation} - - - - - - - - - - - diff --git a/buch/papers/erdbeben/teil1.tex b/buch/papers/erdbeben/teil1.tex index 014b53e..1a893dd 100644 --- a/buch/papers/erdbeben/teil1.tex +++ b/buch/papers/erdbeben/teil1.tex @@ -9,187 +9,213 @@ % (c) 2020 Prof Dr Andreas Müller, Hochschule Rapperswil %% - - \rhead{Kalman-Filter} - \section{Kalman-Filter} -Die interessante Grösse ist also das Integral der Überlagerung zweier Kräfte. -Wir brauchen also die zweite Ableitung der Messung, ohne deren Eigendynamik. -Da wir die äussere Kraft nicht direkt messen können, benötigen wir ein Werkzeug, welches aus der gemessenen Position, die Krafteinwirkung auf unsere System schätzt. -Dies ist eine typische Anwendung für das Kalman-Filter. - -Das Filter schätzt den Zustand eines Systems anhand von Messungen und kann den nächsten Zustand errechnen und aus dieser Schätzung auch eine erwartete Messung herleiten. -Die für das Filter relevante Grösse ist dann nicht mehr die eigentliche Messung, sondern die Differenz aus Messung und Erwartung, da diese Differenz, die Innovation, eine Aussage über die nicht-deterministischen, externen Einflüsse auf das System ermöglicht. -Das Filter berücksichtigt dazu nicht nur die Messung und den Zustand, sondern auch die Unsicherheiten dieser beiden Grössen, welche als Parameter in das Modell des Systems einfliessen. - -Unser Ziel ist es, anhand der Messung die eigentlich interessante Grösse $f$ zu bestimmen. -Dabei wird durch eine deterministische Vorhersage, in dem der Zustand mit der Eigendynamik des Systems multipliziert wird. -Die Idee dahinter ist, dass das Kalman-Filter die nicht-deterministische Grösse $f$ anhand der Messung und der Vorhersage zu bestimmen. - -Für mehrere Dimensionen (x,y,z) würde der Satz von Pythagoras für das System benötigt. -Da sich der Pythagoras bekanntlich nicht linear verhält, kann kein lineares Kalman-Filter implementiert werden. -Da das Kalman-Filter besonders effektiv und einfach für lineare Abläufe geeignet ist, würde eine zweidimensionale Betrachtung den Rahmen dieser Arbeit sprengen. -Einfachheitshalber beschränken wir uns auf den linearen Fall, da dadurch die wesentlichen Punkte bereits aufgezeigt werden. -Für ein nicht-lineares System werden Extended Kalman-Filter benötigt, bei denen die System-Matrix (A) durch die Jacobi-Matrix des System ersetzt wird. - -\subsection{Geschichte} -Das Kalman-Filter wurde 1960 von Rudolf Emil Kalman entdeckt und direkt von der NASA für die Appollo Mission benutzt. -Das Filter kommt mit wenig Rechenleistung aus und war somit dafür geeignet die Rakete bei der Navigation zu unterstützen.
Eine typische Anwendungen des Kalman-Filters ist Glättung von verrauschten Daten und die Schätzung von Parametern. Dies kommt heutzutage in jedem Satellit, Navigationssystem, Smartphones und Videospielen vor. - -\subsection{Wahrscheinlichkeit} -Das Kalman-Filter schätzt den wahrscheinlichsten Wert zwischen Normalverteilungen. -Dies bedeutet, das Filter schätzt nicht nur den Mittelwert, sondern auch die Standartabweichung. -Da Normalverteilungen dadurch vollständig definiert sind, schätzt ein Kalman-Filter die gesamte Verteilungsfunktion des Zustandes. -In der Abbildung~\ref{erdbeben: Zwei Normalverteilungen} sind zwei Funktionen dargestellt. -Die eine Funktion zeigt die errechnete Vorhersage des Zustands, bzw. deren Normalverteilung. -Die andere Funktion zeigt die verrauschte Messung des nächsten Zustand, bzw. deren Normalverteilung. -Wie man am Beispiel der Gauss-Verteilungen in Abblidung~\ref{erdbeben: Zwei Normalverteilungen} sehen kann, ist sowohl der geschätzte Zustand als auch der gemessene Zustand normalverteilt und haben dementsprechend unterschiedliche Standardabweichungen $\sigma$ und Erwartungswerte $\mu$. Dies wird in~\cite{erdbeben:aragher_understanding_2012}beschrieben. +Im letzten Abschnitt haben wir Gleichungen für unser System gefunden. +Als nächstes brauchen wir also ein Werkzeug, +um aus der Messung der Position $s(t)$ den gesammten Zustand $x(t)$ zu schätzen. +Das ist genau das, was Kalman-Filter tun: Ahand von Messungen den Zustand eines Systems schätzen. + +Kalman-Filter wurde 1960 von Rudolf Emil Kalman erfunden und direkt von der NASA für die Appollo Mission benutzt. +Diese Filter kommen mit wenig Rechenleistung aus und waren somit geeignet, die Rakete bei der Navigation zu unterstützen. +Heutige, typische Anwendungen von Kalman-Filtern sind die Glättung verrauschter Daten und die Schätzung von Parametern. +Dies kommt heutzutage in jedem Satellit, Navigationssystem, Smartphones und Videospielen vor. + +Kalman-Filter funktionieren nach folgendem Zwei-Schritt-Verfahren: +Zuerst wird, +ausgehend von der aktuellen Schätzung des Zustands und der Eigendynamik des Systzems, +eine Vorhersage berechnet. +Daraus lässt sich eine erwartete Messung ableiten. +Anschliessend wird diese Vorhersage korrigiert, +wobei die Korrektur abhänging von der Differenz zwischen erwarteter und effektiver Messung ist. + +Dabei sind sowohl die Vorhersage als auch die Messung nur Schätzungen und unweigerlich fehlerbehaftet. +Unter der Annahme, dass die Fehler normalverteilt sind, +lassen sich beide Schätzungen zu einer neuen, optimalen Schätzung kombinieren. +Die genaue Herleitung des Kalman-Filters ist relativ aufwendig +und kann unter Anderem in \cite{erdbeben:skript:wrstat} nachgelesen werden. + +\subsection{Exkurs Wahrscheinlichkeit} +\label{erdbeben:Wahrscheindlichkeit} +Das Kalman-Filter schätzt also den wahrscheinlichsten Wert zwischen zwei Normalverteilungen, +genauer gesagt zwischen einer Messung und einer Vorhersage. +In diesem Abschnitt wollen wir auffrischen, wie dies genau passiert. + +Das Folgende wird in \cite{erdbeben:aragher_understanding_2012} beschrieben. +Wir haben eine Vorhersage aus der Systemdynamik und eine Messung des Zustandes. +Diese widersprechen sich im Allgemeinen. +Jedoch kennen wir auch die Wahrscheinlichkeiten der beiden Aussagen. + \begin{figure} \begin{center} - \includegraphics[width=5cm]{papers/erdbeben/Gausskurve2.pdf} - \caption{Zwei Normalerteilungen; Die eine Funktion zeigt die Vorhersage, die andere die Messung} - \label{erdbeben: Zwei Normalverteilungen} + \includegraphics[width=5cm]{papers/erdbeben/Gausskurve3.pdf} + \caption{ + Seien blau und orange zwei normalverteilte Schätzungen eines Zustandes, etwa eine Vorhersage und eine Messung. + Dann ist die rote Kurve die optimale Schätzung. + Sie entspricht bis auf Normierung dem Produkt von blau und orange.} + \label{erdbeben:Gauss3} \end{center} \end{figure} -Wir haben eine Vorhersage aus der Systemdynamik und eine Messung des Zustandes. -Diese widersprechen sich im Allgemeinen. -Jedoch wissen wir die Wahrscheinlichkeiten der beiden Aussagen. +Abbildung~\ref{erdbeben:Gauss3} zeigt in blau und rot zwei Normalverteilungen, +je eine für die Vorhersage und eine für die Messung. +Diese unterscheiden sich sowohl in ihren Mittelwerten $\mu_{1,2}$, als auch in ihren Standardabweichungen $\sigma_{1,2}$. Um eine genauere Schätzung des Zustandes zu machen, wird nun ein Wert zwischen den beiden Verteilungen berechnet. -Nun wird eine Eigenschaft der Normalverteilung ausgenutzt. Durch das Multiplizieren zweier Normalverteilungen entsteht eine neue Normalverteilung. -Wir haben eine Normalverteilung der Vorhersage: +Nun wird eine Eigenschaft der Normalverteilung ausgenutzt: +Durch das Multiplizieren zweier Normalverteilungen entsteht eine neue Normalverteilung. + +Wir haben also eine Normalverteilung der Vorhersage \[ {y_1}(x;{\mu_1},{\sigma_1})=\frac{1}{\sqrt{2\pi\sigma_1^2}}\quad e^{-\frac{(x-{\mu_1})^2}{2{\sigma_1}^2}} \] -und der Messung: +und der Messung \[ {y_2}(x;{\mu_2},{\sigma_2})=\frac{1}{\sqrt{2\pi\sigma_2^2}}\quad e^{-\frac{(x-{\mu_2})^2}{2{\sigma_2}^2}}. \] -Diesen werden nun multipliziert und durch deren Fläche geteilt um sie wieder zu normieren, $\odot$ beschreibt dabei die Multiplikation und die Normierung auf den Flächeninhalt eins : +Diesen werden nun multipliziert und durch deren Fläche geteilt, +um sie wieder zu normieren. +$\odot$ beschreibt dabei die Multiplikation und die Normierung auf den Flächeninhalt eins: \begin{align*} - {y_f}(x; {\mu_f}, {\sigma_f}) = {y_1}(x;{ \mu_1},{ \sigma_1}) \odot {y_2}(x; {\mu_2}, {\sigma_2}) + {y_f}(x; {\mu_f}, {\sigma_f}) + &= + {y_1}(x;{ \mu_1},{ \sigma_1}) \odot {y_2}(x; {\mu_2}, {\sigma_2}) + \\ &= \frac{1}{\sqrt{2\pi\sigma_1^2}}\quad e^{-\frac{(x-{\mu_1})^2}{2{\sigma_1}^2}} \odot \frac{1}{\sqrt{2\pi\sigma_2^2}}\quad e^{-\frac{(x-{\mu_2})^2}{2{\sigma_2}^2}} \\ &= \frac{ \frac{1}{\sqrt{2\pi\sigma_1^2}}e^{-\frac{(x-{\mu_1})^2}{2{\sigma_1}^2}} \cdot \frac{1}{\sqrt{2\pi\sigma_2^2}}e^{-\frac{(x-{\mu_2})^2}{2{\sigma_2}^2}}}{\int {y_1} {y_2} dx}. \end{align*} -Diese Kombination der beiden Verteilungen resultiert wiederum in einer Normalverteilung -mit Erwartungswert +Die genaue Berechnung ist nicht schwierig aber aufwendig und wird hier deshalb ausgelassen. +Nach einigem Rechnen findet man die Ausdrücke \[ \mu_f = \frac{\mu_1\sigma_2^2 + \mu_2 \sigma_1^2}{\sigma_1^2 + \sigma_2^2} \] -und Varianz +für den neuen Mittelwert und \[ \sigma_f^2 = \frac{\sigma_1^2 \sigma_2^2}{\sigma_1^2 + \sigma_2^2}. \] -Dadurch gleicht sich die neue Kurve den anderen an. Interessant daran ist, dass die fusionierte Kurve sich der genauere Normal-Verteilung anpasst. -Ist ${\sigma_2}$ klein und ${\sigma_1}$ gross, so wird sich die fusionierte Kurve näher an ${y_2}(x;{\mu_2},{\sigma_2})$ begeben. -Somit ist $\mu_f$ ist das gewichtete Mittel der beiden $\mu_{1,2}$, und die Varianzen sind die Gewichte! -Die neue Funktion ist die best mögliche Schätzung für zwei Verteilungen, welche den selben Zustand beschreiben. -Dies ist in der Abbildung~\ref{erdbeben:Gauss3} anhand der rote Funktion ersichtlich. -\begin{figure} - \begin{center} - \includegraphics[width=5cm]{papers/erdbeben/Gausskurve3.pdf} - \caption{Durch das Multiplizieren der blauen und der orangen Verteilung entsteht die die rote, optimale Funktion} - \label{erdbeben:Gauss3} - \end{center} -\end{figure} +für die Varianz. + +Interessant daran ist, dass sich die fusionierte Kurve der genauere Normal-Verteilung anpasst. +Ist ${\sigma_2}$ klein und ${\sigma_1}$ gross, +so wird sich die fusionierte Kurve näher an ${y_2}(x;{\mu_2},{\sigma_2})$ begeben. +$\mu_f$ ist das gewichtete Mittel der beiden $\mu_{1,2}$, und die Varianzen $\sigma_{1,2}$ sind die Gewichte. +Das Interessante an $\mu_{f}$ ist, dass ${\mu_2}$ das Gewicht für ${\sigma_1}$ ist. +Somit ist die Unsicherheit der Messung das Gewicht der Vorhersage und umgekehrt. +Diese neue Funktion ist die best mögliche Schätzung für zwei Verteilungen, welche den selben Zustand beschreiben. +Dies ist in der Abbildung~\ref{erdbeben:Gauss3} anhand der roten Funktion ersichtlich. + Was in zwei Dimensionen erklärt wurde, funktioniert auch in mehreren Dimensionen. Dieses Prinzip mach sich das Kalman Filter zu nutze, und wird von uns für die Erdbeben Berechnung genutzt. -\section{Filter-Matrizen} -Da wir nun ein Werkzeug besitzen, dass die Beschleunigung, welche auf das Gehäuse wirkt, ermitteln kann, wird dieses nun Schritt für Schritt erklärt. -Um den Kalman Filter zu starten, müssen gewisse Bedingungen definiert werden. +\subsection{Filter-Matrizen} +Da wir nun ein Werkzeug besitzen, dass die Beschleunigung, welche auf das Gehäuse wirkt, ermitteln kann, +wird dieses nun Schritt für Schritt erklärt. +Um das Kalman Filter zu starten, müssen gewisse Bedingungen definiert werden. In diesem Abschnitt werden die einzelnen Parameter und Matrizen erklärt und erläutert, wofür sie nützlich sind. -\subsection{Fiter-Agorithmus} -Nachdem alle Parameter aufgestellt sind, wird das Filter initialisiert. -Zuerst wird der nächste Zustand der Masse vorhergesagt, danach wird die Messung präzisiert und laufend aktualisiert. -Das Filter berechnet aufgrund der aktuellen Schätzung eine Vorhersage. -Diese wird, sobald verfügbar, mit der Messung verglichen. -Aus dieser Differenz und den Unsicherheiten des Prozesses ($Q$) und der Messung ($R$) wird der wahrscheinlichste, neue Zustand geschätzt. -Dabei muss genau auf den Index geachtet werden. Nach dem Artikel~\cite{erdbeben:wikipedia} ist die Indexierung so genormt: -Der Zeitschritt wird mit $k$ definiert, $k-1$ ist somit ein Zeitschritt vor $k$. -Auf der linken Seite von | wird der aktuelle Zustand verlangt, bzw. ausgegeben, auf der rechten Seiten den bisherigen Zustand. -Dies bedeutet, dass die Notation $x_{n|m}$ die Schätzung von $x$ zum Zeitpunkt $n$ bis und mit zur Zeitpunkt $m \leq \ n$ präsentiert. +Dabei muss genau auf den Index geachtet werden. +Wir verwenden die Standard-Notation, wie sie auch im Artikel~\cite{erdbeben:wikipedia} zu finden ist. +Sie ist an die Notation der bedingten Wahrscheinlichkeiten angelehnt. +Hierbei steht der betrachtete Zeitschritt links und der gegenwärtige rechts eines Vertikalstrichs. +Dies bedeutet, dass die Notation $x_{n|m}$ die Schätzung von $x$ zum Zeitpunkt $n$ +aufgrund des Wissens bis zum und mit dem Zeitpunkt $m$ repräsentiert. \subsubsection*{Vorhersage} -Im Filterschritt Vorhersage wird der nächste Zustand anhand des Anfangszustand und der Systemmatrix berechnet. -Dies funktioniert mit dem Rechenschritt: +Im Filterschritt Vorhersage wird anhand des aktuellen Zustands und der Systemmatrix eine Schätzung für den nächsten Zustand berechnet. +Die Systemmatrix $A$ aus Gleichung~\eqref{erdbeben:systemmatrix} beschreibt ein kontinuierliches System $\dot x = Ax$. +Wir benötigen jedoch ein Zeit-diskretes System $x_{k+1} = \Phi x_k$. + +Die Exponentialfunktion $\exp(At)$ beschreibt die Entwicklung eine Zustandes im Laufe der Zeit. +Die Übergangs-Matrix $\Phi$ erhalten wir folglich aus der Systemdynamikmatrix $A$ durch die Exponentialfunktion +\[\Phi = \exp(A\Delta t). \] +Die Matrix $\Phi$ beschreibt die Übergänge zwischen zeitlich aufeinanderfolgenden Zuständen $x_{k-1}$ und $x_{k}$ anhand folgender Gleichung: \[ {x_{k|k-1}}=\Phi{x_{k-1|k-1}}= \exp(A\Delta t){x_{k-1|k-1}}. \] -Die Kovarianz $P_{k|k-1}$ wird ebenfalls neu berechnet. Zudem kommt noch die Prozessunsicherheit $Q$ dazu, so dass die Unsicherheit des Anfangsfehlers $P$ laufend verändert. -Dies funktioniert durch multiplizieren der Systemmatrix mit dem aktualisierten Anfangsfehler. -Dazu wird noch die Prozessunsicherheit addiert, somit entsteht die Gleichung +Damit haben wir die Systemdynamik nun in der für unser Kalman-Filter notwendigen Form und können Vorhersagen berechnen. + +Als nächstes benötigen wir die Unsicherheit der Vorhersage. +Im Abschnitt ~\ref{erdbeben:Wahrscheindlichkeit} haben wir dafür die Varianzen der Normalverteilungen verwendet. +Im mehrdimensionalen Fall übernimmt dies die Kovarinanzmatrix $P$. +Sie wird in jedem Schritt aktualisiert. +Hinzu kommt die Prozessunsicherheit $Q$, welche als Parameter in unser Modell einfliesst. +$Q$ beschreibt Unsicherheiten im Modell, +wie etwa unsere Annahme, dass die Kraft sich nicht ändert, +aber auch nicht-modellierbare Einflüsse wie Vibrationen. +$P$ wird dabei laufend aktuallisiert. +Die optimale Gleichung lautet \[ {P_{k|k-1}}=\Phi {P_{k-1|k-1}} {\Phi _{k}}^T + {Q_{k-1}}. \] -Es vergeht genau $\Delta t$ Zeit, und dieser Vorgang wird wiederholt. -Das hochgestellte T bezeichnet die transponierte Matrix. -Dabei wird in den späteren Schritten überprüft, wie genau die letzte Anpassung von $P$ zur Messung stimmt. -Ist der Unterschied klein, wird die Kovarianz $P$ kleiner, ist der Unterschied gross, wird auch die Kovarianz grösser. +Es vergeht genau $\Delta t$ Zeit, und dieser Vorgang wird wiederholt. Das Filter passt sich selber an und korrigiert sich bei grosser Abweichung. \subsubsection*{Messen} -Der Sensor wurde noch nicht benutz, doch genau der liefert Werte für das Filter. -Die aktuellen Messwerte $z$ werden die Innovation $w$ mit dem Zustandsvektor $x$ und der Messmatrix $H$ zusammengerechnet. -Hier bei wird lediglich die Messung mit dem Fehler behaftet, und die Messmatrix $H$ mit der Vorhersage multipliziert. +Der Sensor wurde noch nicht benutz, doch genau der liefert die Messwerte $z_k$ für unser Filter. +Aus der Vorhersage des Zustandes $x_{k|k-1}$ und der Messmatrix $H$ erhalten wird eine Vorhersage der Messung. +Die Innovation \[ -{w_{k}}={z_{k}}-{H}{x_{k|k-1}}. +{w_{k}}={z_{k}}-{H}{x_{k|k-1}} \] -Die Innovation ist der Teil der Messung, die nicht durch die Systemdynamik erklärt werden kann. -Die Hilfsgröße Innovation beschreibt, wie genau die Vorhersage den aktuellen Messwert mittels der Systemmatrix $\Phi$ beschreiben kann. +beschreibt, wie genau die Vorhersage $x_{k|k-1}$ zur aktuellen Messung $z_k$ passt. +Die Innovation ist also derjenige Teil der Messung, der nicht im Modell erfasst ist. +Dies leuchtet ein, eine Innovation von $0$ bedeutet, dass die Messung nichts Neues hervorbrachte. Für eine schlechte Vorhersage wird die dazugehörige Innovation gross, für eine genaue Vorhersage dagegen klein sein. -Entsprechende Korrekturen müssen dann gross bzw. nur gering ausfallen. -Innovation = Messung - Vorhersage. Dies leuchtet ein, eine Innovation von 0 bedeutet, dass die Messung nichts Neues hervorbrachte. - -Im nächsten Schritt wir analysiert, mit welcher Kovarianz weiter gerechnet wird. -Hierbei wird die Unsicherheit $P$, die Messmatrix $H$ und die Messunsicherheit $R$ miteinander verrechnet. -\[ -{S_{k}}={H}{P_{k|k-1}}{H}^T+{R_{k}} -\] +Entsprechende Korrekturen werden dann gross bzw. nur gering ausfallen. \subsubsection*{Aktualisieren} -Im nächsten Schritt kommt nun die Wahrscheinlichkeit dazu. -\[{K_{k}}= {P_{k|k-1}} {H^T}{S_{k}^{-1}}\] -Die Grösse $K$ wird Kalman-Gain genannt. -Das Kalman-Gain gibt dem Zustand die Gewichtung, bzw. wie die Vorhersage auf den Zustand passt. -Vereinfacht gesagt: Es wird das das Verhältnis zwischen der Unsicherheit der Vorhersage $P_k$ zu der zugehörigen Messunsicherheit $R_k$ gebildet. -In unserem Fall wird werden die Elemente der Kalman-Matrix vorweg berechnet, da das Kalman-Gain ohne Messungen auskommt. - -Anhand der Informationen aus der Innovation wird das Kalman-Gain $K$ gebildet. Dabei beschreibt das Kalman-Gain die Wirkung der Innovation auf den geschätzten Zustand. So wird das System aktualisiert. + +Für eine optimale Schätzung des Zustandes muss die Vorhersage entsprechend der Innovation korrigiert werden. +In der Literatur findet man für eine optimales Korrektur die Gleichungen +\begin{align*} +{S_{k}} &={H}{P_{k|k-1}}{H}^T+{R_{k}} +\\ +{K_{k}} &= {P_{k|k-1}} {H^T}{S_{k}^{-1}} +\end{align*} +Dabei ist $K$ das Kalman-Gain. +$K$ beschreibt, wie die Vorhersage korrigiert werden muss. +Die optimale Schätzung des neuen Zustandes wird dann zu \[ -{x_{k|k}}={x_{k|k-1}}+{K_{k}}{w_{k}} +{x_{k|k}}={x_{k|k-1}}+{K_{k}}{w_{k}}. \] -Dabei wird der Unterschied zwischen dem erwarteten, errechneten, Zustand und dem gemessenen Zustand berechnet. - -Dazu kommt eine neue Kovarianz für den nächste Vorhersageschritt: +Dazu kommt eine neue Kovarianz $P$ für den nächste Vorhersageschritt: \[ -{P_{k|k}}=(I-{K_{k}}{H}){P_{k|k-1}} +{P_{k|k}}=(I-{K_{k}}{H}){P_{k|k-1}}. \] -Der ganze Algorithmus und beginnt wieder mit der Vorhersage +Der ganze Algorithmus ist nun vollständig und beginnt wieder mit der Vorhersage \[ -{x_{k|k-1}}=\Phi{x_{k-1|k-1}}= \exp(A\Delta t){x_{k|k-1}}. +{x_{k+1|k}}=\Phi{x_{k|k}}= \exp(A\Delta t){x_{k|k}}. \] -\subsection{Anfangsbedingungen} +\subsection{Parameter und Anfangsbedingungen} +Die Grössen $P$, $Q$, $R$ und $\Phi$ können grundsätzlich in jedem Zeitschritt ändern. +Für die meisten Anwendungen sind sie jedoch konstant und fliessen als Parameter ins Modell ein. +Aufgrund der iterativen Arbeitsweise von Kalman-Filtern benötigen wir zudem ein paar Anfangswerte. + \subsubsection*{Anfangszustand $x$} -Das Filter benötigt eine Anfangsbedingung. +Für die erste Vorhersage benötigt das Filter einen Anfangszustand. In unserem Fall ist es die Ruhelage, die Masse bewegt sich nicht. Zudem erfährt die Apparatur keine äussere Kraft. \[ {x_0 }= \left( \begin{array}{c} {s_0}\\ {v_0}\\{f_0}\end{array}\right) = \left( \begin{array}{c} 0\\ 0\\ 0\end{array}\right) \] +\subsubsection*{Systemmatrix $A$ und $\Phi$} +Für unseren Seismographen haben wir die entsprechende Matrixdarstellung +in Gleichung~\eqref{erdbeben:systemmatrix} bereits gefunden. +Zudem haben wir weiter oben bereits beschrieben, +wie wir mittels Exponentialfunktion zu einer zeitdiskreten Beschreibung für das Kalman-Filter kommen. +Es gilt +\[ \Phi = \exp(A \Delta t) .\] + \subsubsection*{Anfangsfehler / Kovarianzmatrix $P$} -Da auch der Anfangszustand fehlerhaft sein kann, wird für das Filter ein Anfangsfehler verwendet. +Da auch der Anfangszustand fehlerhaft sein kann, wird für das Filter eine Anfangsunsicherheit verwendet. Auf der Diagonalen werden die Varianzen eingesetzt, in den restlichen Felder stehen die Kovarianzen. -Zur Erinnerung: Die Varianz ist ein Mass für die Streuung eines Wertes, die Kovarianz hingegen beschreibt die Abhängigkeit der Streuungen zweier Werte. - -Kovarianz: Cov(x, y) und Varianz: Var(x) = Cov(x, x) +Für einen gut bekannten Zustandsvektor können kleine Werte eingesetzt werden, für ungenaue Anfangsbedingungen sollten grosse Werte verwendet werden. +Grosse Werte ermöglichen dem Filter sich schnell einzupendeln. In unserem Fall ist der Anfangszustand gut bekannt. -Wir gehen davon aus, dass das System in Ruhe und in Abwesenheit eines Erdbeben startet, somit kann die Matrix mit Nullen bestückt werden. -Als Initialwert für die Kovarianzmatrix ergibt sich +Wir gehen davon aus, +dass das System in Ruhe und in Abwesenheit eines Erdbeben startet. +Somit kann die Matrix mit Nullen bestückt werden und wir starten mit \[ {P_0 }= \left( @@ -200,35 +226,13 @@ Als Initialwert für die Kovarianzmatrix ergibt sich \end{array} \right). \] -Diese Matrix beschreibt die Unsicherheit des geschätzten Zustandes und wird sowohl für die Vorhersage als auch die Korrektur benötigt. -Sie wird nach jeder Schätzung aktualisiert. -Für einen gut bekannten Zustandsvektor können kleine Werte eingesetzt werden, für ungenaue Anfangsbedingungen sollten grosse Werte verwendet werden. -Grosse Werte ermöglichen dem Filter sich schnell einzupendeln. -\subsubsection*{Dynamikmatrix $A$} -Das Kalman-Filter benötigt für die Vorhersage des nächsten Zustandes eine Beschreibung der Systemdynamik. -Die Dynamikmatrix bildet den Kern des Filters. Diese wurde weiter oben bereits beschrieben. -Dabei wollen wird die äussere Kraft des Systems ermitteln. -Da nichts über die äussere Kraft bekannt ist, müssen wir annehmen das deren Ableitung 0 ist. -Die System-Matrix lautet daher: -\[ -A = \left( - \begin{array}{ccc} -0 & 1& 0 \\ -- \frac{D}{m} &-\frac{2k}{m} & \frac{1} {m}\\ -0 & 0& 0\\ -\end{array}\right) - \] -Dabei soll der Kalman-Filter in diskreten Zeitschritten $\Delta t$ arbeiten. -$A$ beschreibt ein kontinuierliches System ($\dot x = Ax$), wir benötigen jedoch ein Zeit-diskretes System $x_{k+1} = \Phi x_k$. -Die Übergangs-Matrix erhalten wir aus der Systemdynamikmatrix mittels Exponentialfunktion: -\[\Phi = \exp(A\Delta t). \] -Die Matrix $\Phi$ beschreibt die Übergänge zwischen zeitlich aufeinanderfolgenden Zuständen $x_{k-1}$ und $x_{k}$ \subsubsection*{Prozessrauschkovarianzmatrix $Q$} Die Prozessrauschmatrix teilt dem Filter mit, wie sich der Prozess verändert. Die Matrix $Q$ beschreibt die Unsicherheit, die der Prozess mit sich bringt. -Bei unserem Modell könnte das beispielsweise ein Windstoss an die Masse sein oder auch die Ungenauigkeiten im Modell, wie die Annahme das dich die Kraft nicht ändert. +Bei unserem Modell könnte das beispielsweise ein Windstoss an die Masse sein +oder auch Ungenauigkeiten im Modell, wie die Annahme, dass sich die Kraft nicht ändert. Für uns wäre dies: \[ Q = \left( @@ -241,9 +245,9 @@ Q = \left( Die Standabweichungen müssten statistisch ermittelt werden, da der Fehler nicht vom Sensor kommt und somit nicht vom Hersteller gegeben ist. \subsubsection*{Messmatrix $H$} -Die Messmatrix gibt an, welche Parameter gemessen werden. -$H$ ist die Matrix für die Vorhersage der Messung. -In unserem Falle ist es die Position der Massen. +Die Messmatrix gibt an, welche Zustände gemessen werden. +$H$ ist die Matrix, welche aus der Vorhersage des Zustand eine Vorhersage der Messung erzeugt. +In unserem Falle messen wir nur die Position der Massen und verwenden deshalb \[ H = (1, 0, 0) \] @@ -252,49 +256,48 @@ H = (1, 0, 0) Die Messrauschkovarianzmatrix beinhaltet, wie der Name schon sagt, das Rauschen der Messung. In unserem Fall wird nur die Position der Masse gemessen. Da wir keine anderen Sensoren haben ist $R$ lediglich: \[ -R= ({\sigma_\mathrm{sensor}}^2). +R= (\sigma_\mathrm{sensor}^2). \] Diese Messrauchen wird meistens vom Sensorhersteller angegeben. -Für unsere theoretische Apparatur wird hier ein kleiner Fehler eingesetzt da heutige Sensoren sehr genau messen können. +Für unsere theoretische Apparatur wird hier ein kleiner Fehler eingesetzt, +da heutige Sensoren sehr genau messen können. +\clearpage \subsection{Zusammenfassung } -Zusammenfassend kann das Kalman-Filter in offizieller Typus dargestellt werden. -Dabei beginnt das Filter mit dem Anfangszustand für $k=0$ - -1. Nächster Zustand vorhersagen +Das Filter beginnt mit dem Anfangszustand für $k=0$. +Anschliessend werden folgende Schritte iterativ ausgeführt: +\begin{enumerate} +\item Nächster Zustand vorhersagen \[ -{x_{k|k-1}}=\Phi{x_{k-1|k-1}}= \exp(A\Delta t){x_{k-1|k-1}}. +{x_{k|k-1}}=\Phi{x_{k-1|k-1}}= \exp(A\Delta t){x_{k-1|k-1}} \] -2. Nächste Fehlerkovarianz vorhersagen + \item Nächste Fehlerkovarianz vorhersagen \[ -{P_{k|k-1}}=\Phi {P_{k-1|k-1}} {\Phi _{k}}^T + {Q_{k-1}}. +{P_{k|k-1}}=\Phi {P_{k-1|k-1}} {\Phi _{k}}^T + {Q_{k-1}} \] -3. Zustand wird gemessen +\item Innovation (= Messung - Vorhersage) \[ -{w_{k}}={z_{k}}-{H}{x_{k|k-1}}. +{w_{k}}={z_{k}}-{H}{x_{k|k-1}} \] -4. Innovation (= Messung - Vorhersage) -\[ -{S_{k}}={H}{P_{k|k-1}}{H}^T+{R_{k}} -\] - -5. Das Kalman Filter anwenden -\[ -{K_{k}}= {P_{k|k-1}} {H^T}{S_{k}^{-1}} -\] +\item Optimales Kalman-Gain berechnen +\begin{align*} +{S_{k}} &={H}{P_{k|k-1}}{H}^T+{R_{k}}\\ +{K_{k}} &= {P_{k|k-1}} {H^T}{S_{k}^{-1}} +\end{align*} -6. Schätzung aktualisieren +\item Schätzung aktualisieren \[ {x_{k|k}}={x_{k|k-1}}+{K_{k}}{w_{k}} \] -7. Fehlerkovarianz aktualisieren +\item Fehlerkovarianz aktualisieren \[ {P_{k|k}}=(I-{K_{k}}{H}){P_{k|k-1}} \] -8. Die Outputs von $k$ werden die Inputs für ${k-1}$ und werden wieder im Schritt 1 verwendet +\end{enumerate} +Die Outputs von $k$ werden die Inputs für ${k+1}$ und werden wieder in Schritt 1 verwendet. diff --git a/buch/papers/munkres/teil1.tex b/buch/papers/munkres/teil1.tex index 0c23c7c..5f2fed6 100644 --- a/buch/papers/munkres/teil1.tex +++ b/buch/papers/munkres/teil1.tex @@ -18,6 +18,8 @@ Um dieses Problem in einer einfachen, händischen Art und Weise zu lösen wurde Als Beispiel betrachten wir den Fall, wo ein Bauunternehmer einen Bauingenieur beauftragt, eine optimale Transportroute für die Umplatzierung seiner Kräne zu eruieren. Das heisst, die Transportstrecke für die Umplatzierung seine Kräne soll möglichst klein werden. Die Frage lautet: Wie sind die Kräne umzusetzen, damit deren Transportstrecke minimal wird? Bei der normalen Optimierung dürfen normalerweise beliebige reelle Werte $\mathbb{R}$ angenommen werden. +\index{Kran}% +\index{Baustelle}% \index{Optimierung}% Für solche Optimierungsprobleme für reelle Variablen sind verschiedene Verfahren entwickelt worden, die im Allgemeinen auch sehr effizient sind. Das reelle Problem ist also in einer einfachen Art und Weise lösbar. |