09_interpolation_unregelmaessig_verteilter_daten.tex

\chapter{%
    Interpolation unregelmäßig verteilter Daten%
}

\section{%
    \name{Voronoi}-Diagramm und \name{Delaunay}-Triangulierung%
}

\textbf{\name{Voronoi}-Diagramm}:
Gegeben sei eine Menge $\S := \{\vecs{x}{1}, \dotsc, \vecs{x}{N}\} \subset \real^n$ von
$N$ Punkten und eine Metrik $d$ auf $\real^n$.
Die \begriff{\name{Voronoi}-Zelle} von $\vecs{x}{i}$ ist definiert als\\
$V(\vecs{x}{i}) := \{\vec{x} \in \real^n \;|\;
\forall_{j\not=i}\; d(\vec{x}, \vecs{x}{i}) < d(\vec{x}, \vecs{x}{j})\}$.\\
Das \begriff{\name{Voronoi}-Diagramm} von $\S$ ist $\Vor(\S) := \{V(\vecs{x}{i})\}_{i=1}^N$.

\textbf{\name{Delaunay}-Triangulierung}:
Die \begriff{\name{Delaunay}-Triangulierung} $\Del(\S)$ ist dual zu $\Vor(\S)$.\\
Zwei Punkte $\vecs{x}{i}$, $\vecs{x}{j}$ werden verbunden, falls
$V(\vecs{x}{i})$ und $V(\vecs{x}{j})$ eine gemeinsame Kante besitzen.
$\vecs{x}{i}, \vecs{x}{j}, \vecs{x}{k}$ bilden ein Dreieck genau dann, wenn
kein anderer Punkt aus $\S$ im Umkreis des Dreiecks
$\Delta\vecs{x}{i}\vecs{x}{j}\vecs{x}{k}$ liegt.
$\Del(\S)$ maximiert lexikografisch die minimalen Innenwinkel der Triang. und außerdem
das Verhältnis $a := \frac{r_\text{Inkreis}}{r_\text{Umkreis}}$.

\textbf{\name{Bowyer}-\name{Watson}-Algorithmus}:
Der \begriff{\name{Bowyer}-\name{Watson}-Algorithmus} konstruiert $\Del(\S)$ inkrementell.
Um einen Punkt $\vecs{x}{i}$ in die aktuelle Delaunay-Triang. der Punkte
$\vecs{x}{1}, \dotsc, \vecs{x}{i-1}$ einzufügen,
wird zunächst betrachtet, in welchen Dreiecksumkreisen $\vecs{x}{i}$ liegt.
Diese Dreiecke werden dann aus der aktuellen Delaunay-Triang. gelöscht,
was ein sternförmiges Gebiet erzeugt.
Verbinde nun den Punkt $\vecs{x}{i}$ mit allen Ecken, um die neue Delaunay-Triang. zu erhalten.

\linie

Seien nun Datenwerte $f_1, \dotsc, f_N$ an den $\vecs{x}{1}, \dotsc, \vecs{x}{N}$ gegeben.

\textbf{Interpolation mit \name{Voronoi}-Zerlegung}:
Mit der Voronoi-Zerlegung kann man interpolieren,
indem man jeder Zelle $V(\vecs{x}{i})$ den Wert $f_i$ zuweist
(nächster Nachbar).

\textbf{Interpolation mit \name{Delaunay}-Triangulierung}:
Mit der Delaunay-Triang. kann man interpolieren,
indem man lineare Interpolation in jedem Dreieck durchführt.

\section{%
    \name{Shepard}-Interpolation%
}

\textbf{\name{Shepard}-Interpolation}:
Die \begriff{\name{Shepard}-Interpolation} (oder \begriff{inverse Distanzwichtung})
interpoliert $f_i \in \real$, $\vecs{x}{i} := (x_i, y_i) \in \real^2$ ($i = 1, \dotsc, N$)
mit dem Ansatz $f(x, y) := \sum_{j=1}^N f_j w_j(x, y)$.
Damit erhält man $\forall_{i=1,\dotsc,N}\; f_i = f(x_i, y_i) \iff
\forall_{i,j=1,\dotsc,N}\; w_j(x_i, y_i) = \delta_{j,i}$.\\
Dies motiviert den Ansatz
$w_j(x, y) := \frac{\varphi_j(x, y)^\mu}{\sum_{k=1}^N \varphi_k(x, y)^\mu}$ mit
$\varphi_j(x, y) := \frac{1}{|\vec{x} - \vecs{x}{j}|}$ und $\mu \in (0, \infty)$.\\
Es gilt $w_j(x, y) \in [0, 1]$, $\sum_{j=1}^N w_j(x, y) = 1$ und $w_j(x_i, y_i) = \delta_{j,i}$\\
(wenn man $\varphi_j(x_j, y_j) := \infty$ und $w_j(x_j, y_j) := 1$ setzt).

\linie

\textbf{lokale \name{Shepard}-Interpolation}:
Bei der \begriff{lokalen \name{Shepard}-Interpolation} verwendet man\\
Gewichtsfunktionen $w_j$, die nur von einer lokalen Teilmenge von Samplepunkten abhängen.
Definiert man $r_j = r_j(x, y) := |\vec{x} - \vecs{x}{j}|$, dann
kann man z.\,B.
\begin{itemize}
    \item
    $\varphi_j(x, y) := \left(\frac{R - r_j}{R \cdot r_j}\right)^2$ für $r_j \in [0, R]$
    (\begriff{modifizierte \name{Shepard}-Gewichte}),

    \item
    $\varphi_j(x, y) := \frac{1}{r_j}$ für $r_j \in [0, \frac{R}{3}]$ und
    $\varphi_j(x, y) := \frac{27}{4R} \left(\frac{r_j}{R} - 1\right)^2$ für
    $r_j \in [\frac{R}{3}, R]$ oder alternativ

    \item
    $\varphi_j(x, y) := 1 - \frac{r_j}{R}$ für $r_j \in [0, R]$
    (\begriff{\name{Franke}-\name{Little}-Gewichte}) verwenden.
\end{itemize}

\pagebreak

\section{%
    Methode der radialen Basisfunktionen%
}

\textbf{Methode der radialen Basisfunktionen}:
Seien $N$ paarw. disj. Datenpunkte $\vecs{x}{i} \in \real^n$ und Werte $f_i \in \real$
für $i = 1, \dotsc, N$ gegeben.
Die \begriff{Methode der radialen Basisfunktionen} arbeitet ähnlich wie die Shepard-Intp.,
mit dem Unterschied, dass hier $f(\vec{x}) := \sum_{i=1}^N \lambda_i \phi(|\vec{x} - \vecs{x}{i}|)$
mit einer \begriff{radialen Basisfunktion (RBF)} $\phi(r)\colon [0, \infty) \to \real$
angesetzt wird (d.\,h. die Koef"|fizienten sind jetzt allgemein und
die Funktionen sind zwingend radialsymmetrisch).

\textbf{gebräuchliche Basisfunktionen}:
$\phi(r) := e^{-(\varepsilon r)^2}$
(\begriff{\name{Gauß} (GA)}),\\
$\phi(r) := \frac{1}{1 + (\varepsilon r)^2}$
(\begriff{invers quadratisch (IQ)}),
$\phi(r) := \frac{1}{\sqrt{1 + (\varepsilon r)^2}}$
(\begriff{invers multiquadratisch (IMQ)}),\\
$\phi(r) := \sqrt{1 + (\varepsilon r)^2}$
(\begriff{\name{Hardy}-multiquadratisch (MQ)}),
$\phi(r) := r$
(\begriff{linear}),
$\phi(r) := r^3$
(\begriff{kubisch}),\\
$\phi(r) := r^2 \ln r$
(\begriff{Thin Plate Spline (TPS)})
% \begin{itemize}
%     \item
%     $\phi(r) := e^{-(\varepsilon r)^2}$
%     (\begriff{\name{Gauß} (GA)}),
%
%     \item
%     $\phi(r) := \frac{1}{1 + (\varepsilon r)^2}$
%     (\begriff{invers quadratisch (IQ)}),
%
%     \item
%     $\phi(r) := \frac{1}{\sqrt{1 + (\varepsilon r)^2}}$
%     (\begriff{invers multiquadratisch (IMQ)}),
%
%     \item
%     $\phi(r) := \sqrt{1 + (\varepsilon r)^2}$
%     (\begriff{\name{Hardy}-multiquadratisch (MQ)}),
%
%     \item
%     $\phi(r) := r$
%     (\begriff{linear}),
%
%     \item
%     $\phi(r) := r^3$
%     (\begriff{kubisch}),
%
%     \item
%     $\phi(r) := r^2 \ln r$
%     (\begriff{Thin Plate Spline (TPS)}).
% \end{itemize}

\textbf{Ermitteln der Koef"|fizienten}:
Die Koef"|fizienten bekommt man aus dem LGS\\
$\forall_{i=1,\dotsc,N}\; f(\vecs{x}{i}) = f_i
\iff \Phi\vec{\lambda} = \vec{f}$
mit $\Phi := (\phi_{i,j})_{i,j=1}^N$, $\phi_{i,j} := \phi(|\vecs{x}{i} - \vecs{x}{j}|)$,
$\vec{\lambda} := (\lambda_j)_{j=1}^N$ und $\vec{f} := (f_i)_{i=1}^N$.
Ist $\Phi$ invertierbar, so gilt $\vec{\lambda} = \Phi^{-1} \vec{f}$.

\linie

Im Folgenden werden hinreichende Bedingungen für die Invertierbarkeit von $\Phi$
angegeben.

\textbf{vollständig monoton}:
Eine Funktion $\psi\colon [0, \infty) \to \real$ heißt \begriff{vollständig monoton}, falls
\begin{itemize}
    \item
    $\psi \in \C^0([0, \infty))$,

    \item
    $\psi \in \C^\infty((0, \infty))$ und

    \item
    $\forall_{\ell \in \natural_0} \forall_{r > 0}\; (-1)^\ell \psi^{(\ell)}(r) \ge 0$.
\end{itemize}

\textbf{Satz von \name{Schoenberg}}:
Sei $\psi(r) := \phi(\sqrt{r})$.
Ist $\psi$ vollständig monoton, aber nicht konstant auf $[0, \infty)$,
dann ist $\Phi := (\phi_{i,j})_{i,j=1}^N$ für paarw. disjunkte Punkte
$\vecs{x}{i}$ p.\,d. und insb. invertierbar.

\textbf{Beispiel}:
Für die Gauß-RBF ist $\psi(r) = e^{-\varepsilon^2 r}$
und Schoenberg ist anwendbar
($(-1)^\ell \psi^{(\ell)}(r) > 0$).\\
Für die MQ-RBF ist $\psi(r) = \sqrt{1 + \varepsilon^2 r}$
und Schoenberg ist nicht anwendbar
($-\psi'(r) < 0$).

\linie

\textbf{Satz 1 von \name{Micchelli}}:
Sei $\psi(r) := \phi(\sqrt{r})$.
Ist $\psi \in \C^0([0, \infty))$,
$\forall_{r>0}\; \psi(r) > 0$ und
$\psi'$ vollständig monoton, aber nicht konstant auf $[0, \infty)$,
dann ist $\Phi := (\phi_{i,j})_{i,j=1}^N$ für paarw. disjunkte Punkte
$\vecs{x}{i}$ invertierbar.

\textbf{Beispiel}:
Für die MQ-RBF ist der Satz anwendbar.
Allerdings kann man den Satz für TPS nicht anwenden,
weil für $\psi(r) = r \ln\sqrt{r}$ gilt,
dass $\psi(r) < 0$ für $r > 0$ klein
(Schoenberg geht auch nicht, weil $-\psi'(r) = -\frac{1}{2} - \ln\sqrt{r} < 0$
für $r$ groß).

\linie

\textbf{augmentierte RBF-Interpolation}:
Sei $\{p_k\}_{k=1}^M$ eine Basis von $\PP_m(\real^d)$
($d$-variate Polynome vom Grad $\le m$, $M := \binom{m+d}{m}$) und $m \in \natural_0$.
Dann ist der Ansatz für die \begriff{augmentierte RBF-Interpolation}
$f(\vec{x}) = \sum_{i=1}^N \lambda_j \phi(|\vec{x} - \vecs{x}{j}|) +
\sum_{k=1}^M \gamma_k p_k(\vec{x})$.
Damit die Interpolation nicht unterbestimmt ist, legt man die
zusätzlichen Bedingungen $\forall_{k=1,\dotsc,M}\; \sum_{j=1}^N \lambda_j p_k(\vecs{x}{j}) = 0$
fest.\\
Man erhält das Interpolations-LGS
$A\smallpmatrix{\vec{\lambda}\\\vec{\gamma}} = \smallpmatrix{\vec{f}\\\vec{0}}$ mit
$A := \smallpmatrix{\Phi&P\\P^\tp&0}$, $\Phi := (\phi(|\vecs{x}{j} - \vecs{x}{k}|))_{j,k=1}^N$ und
$P := (p_k(\vecs{x}{j}))_{j,k=1}^{N,M}$.

\textbf{Satz 2 von \name{Micchelli}}:
Seien $\psi(r) := \phi(\sqrt{r})$ und $m \in \natural_0$.
Ist $\psi \in \C^0([0, \infty))$,
$\psi^{(m+1)}$ vollständig monoton, aber nicht konstant auf $[0, \infty)$, und
$\Rang(P) = M$ (d.\,h. $P$ hat vollen Spaltenrang),
dann ist $A$ für paarw. disjunkte Punkte $\vecs{x}{i}$ invertierbar.

\pagebreak