Remez-Algorithmus

Der Remez-Algorithmus in der Approximationstheorie ist ein Minimax-Approximations-Algorithmus. Als solcher minimiert er die maximale absolute Differenz zwischen dem gesuchten Polynom p {\displaystyle p} (vorgegebenen Maximalgrades n {\displaystyle n} ) und der gegebenen (in einem Intervall) stetigen Funktion f {\displaystyle f} . Er ist benannt nach dem sowjetischen Mathematiker Jewgeni Jakowlewitsch Remes, der ihn im Jahr 1934 vorgestellt hat.^[1]

Der Algorithmus setzt dabei ganz wesentlich auf den Alternantensatz von Pafnuti Lwowitsch Tschebyschow, indem er iterativ die genannte Differenz an n + 2 {\displaystyle n+2} Stützstellen im Intervall auswertet und daraus neue n + 2 {\displaystyle n+2} Stützstellen berechnet.

Gegeben

Ein Intervall [ a , b ] ⊂ R {\displaystyle [a,b]\subset \mathbb {R} } und eine stetige Funktion f : [ a , b ] → R {\displaystyle f\colon [a,b]\to \mathbb {R} } ; ferner ein maximaler Polynomgrad n {\displaystyle n} .

Gesucht

Ein Polynom p ( x ) ∈ R [ x ] {\displaystyle p(x)\in \mathbb {R} [x]} vom Grad höchstens n {\displaystyle n} , welches

max a ≤ x ≤ b | f ( x ) − p ( x ) | {\displaystyle \max _{a\leq x\leq b}|f(x)-p(x)|}

minimiert.

Aus dem Alternantensatz folgt, dass das Polynom im Limes eindeutig bestimmt ist und dass es n + 2 {\displaystyle n+2} Punkte ( a ≤ ) x 0 < x 1 < … < x n + 1 ( ≤ b ) {\displaystyle (a\leq )\;x_{0}<x_{1}<\ldots <x_{n+1}\;(\leq b)} gibt, für die gilt

f ( x i ) − p ( x i ) = ± ( − 1 ) i ⋅ E {\displaystyle f(x_{i})-p(x_{i})=\pm (-1)^{i}\cdot E}

mit der Abweichung E := max a ≤ x ≤ b | f ( x ) − p ( x ) | = ‖ f − p ‖ ∞ {\displaystyle E:=\max _{a\leq x\leq b}|f(x)-p(x)|=\|f-p\|_{\infty }} (Supremumsnorm).

Das gesuchte Polynom sei mit

p ( x ) =: c 0 + c 1 ⋅ x + ⋯ + c n ⋅ x n {\displaystyle p(x)=:c_{0}+c_{1}\cdot x+\,\cdots \,+c_{n}\cdot x^{n}}

bezeichnet. Es gilt also, die n + 1 {\displaystyle n+1} Unbekannten

c 0 , c 1 , … , c n {\displaystyle c_{0},c_{1},\ldots ,c_{n}}

zu bestimmen.

Man beginnt mit den Extremstellen des Tschebyschow-Polynoms erster Art vom Grad n + 1 {\displaystyle n+1}

x i = a + b 2 − b − a 2 ⋅ cos ⁡ ( i π n + 1 ) {\displaystyle x_{i}={\frac {a+b}{2}}-{\frac {b-a}{2}}\cdot \cos \left({\frac {i\pi }{n+1}}\right)}   mit   i = 0 , … , n + 1 {\displaystyle i=0,\ldots ,n+1} .

Ein solcher Satz von Punkten wird „Referenz“^[2] genannt. Es ist

a = x 0 < x 1 < . . . < x n + 1 = b {\displaystyle a=x_{0}<x_{1}<\,...\,<x_{n+1}=b} .

Gesucht wird das Polynom p ( x ) {\displaystyle p(x)} vom Grad ≤ n {\displaystyle \leq n} , dessen Fehlerfunktion f ( x ) − p ( x ) {\displaystyle f(x)-p(x)} auf der im vorangegangenen Schritt gewonnenen Referenz x 0 < x 1 < ⋯ < x n + 1 {\displaystyle x_{0}<x_{1}<\dots <x_{n+1}} alterniert. D. h. gesucht sind

c 0 , c 1 , … , c n {\displaystyle c_{0},c_{1},\ldots ,c_{n}}   und   E {\displaystyle E} .

mit

p ( x i ) + ( − 1 ) i ⋅ E = f ( x i ) {\displaystyle p(x_{i})+(-1)^{i}\cdot E=f(x_{i})}   für   i = 0 , … , n + 1 {\displaystyle i=0,\ldots ,n+1} .

Diese Aufgabe hat als Eingabe nur die Referenz und von der Funktion f {\displaystyle f} nur die n + 2 {\displaystyle n+2} Werte f ( x i ) {\displaystyle f(x_{i})} . Sie kann auch als Optimierungsaufgabe formuliert werden, nämlich als beste Approximation auf der Referenz, und ist eindeutig lösbar.

Das zu lösende Gleichungssystem in Matrixdarstellung:^[3]

( 1 x 0 ⋯ x 0 n 1 1 x 1 ⋯ x 1 n − 1 ⋮ ⋮ ⋮ ⋮ ⋮ 1 x n ⋯ x n n ( − 1 ) n 1 x n + 1 ⋯ x n + 1 n ( − 1 ) n + 1 ) ( c 0 c 1 ⋮ c n E ) = ( f ( x 0 ) f ( x 1 ) ⋮ f ( x n ) f ( x n + 1 ) ) {\displaystyle \left({\begin{matrix}1&x_{0}&\cdots &x_{0}^{n}&1\\1&x_{1}&\cdots &x_{1}^{n}&-1\\\vdots &\vdots &\vdots &\vdots &\vdots \\1&x_{n}&\cdots &x_{n}^{n}&(-1)^{n}\\1&x_{n+1}&\cdots &x_{n+1}^{n}&(-1)^{n+1}\\\end{matrix}}\right){\begin{pmatrix}c_{0}\\c_{1}\\\vdots \\c_{n}\\E\end{pmatrix}}={\begin{pmatrix}f(x_{0})\\f(x_{1})\\\vdots \\f(x_{n})\\f(x_{n+1})\\\end{pmatrix}}}

Damit hat man n + 1 {\displaystyle n+1} neue Koeffizienten

c 0 , c 1 , … , c n {\displaystyle c_{0},c_{1},\ldots ,c_{n}}

für das Polynom p ( x ) {\displaystyle p(x)} und eine neue »Referenzabweichung«  E {\displaystyle E} .

Ausgehend vom Polynom p ( x ) {\displaystyle p(x)} gilt es, n + 2 {\displaystyle n+2} Extremstellen x ~ 0 , x ~ 1 , … , x ~ n + 1 {\displaystyle {\tilde {x}}_{0},{\tilde {x}}_{1},\ldots ,{\tilde {x}}_{n+1}} der Fehlerfunktion

f ( x ) − p ( x ) {\displaystyle f(x)-p(x)}

zu finden. Ist f {\displaystyle f} differenzierbar, kann man Extremstellen oft durch Nullsetzen der Ableitungsfunktion f ′ − p ′ {\displaystyle f'-p'} gewinnen.

In jedem Fall lässt sich das Intervall [ a , b ] {\displaystyle [a,b]} in die durch die aktuelle Fehlerfunktion spezifizierbaren Teilintervalle [ s i , s i + 1 ] {\displaystyle [s_{i},s_{i+1}]} folgendermaßen aufteilen. Da mit f {\displaystyle f} auch f − p {\displaystyle f-p} stetig ist, gibt es nach dem Zwischenwertsatz für alle i = 1 , … , n + 1 {\displaystyle i=1,\ldots ,n+1} Nullstellen s i {\displaystyle s_{i}} der Fehlerfunktion (in der Abbildung Schnittstellen der roten Funktion f {\displaystyle f} mit dem blauen Polynom p {\displaystyle p} ) im Intervall [ x i − 1 , x i ] {\displaystyle [x_{i-1},x_{i}]} , da das Vorzeichen der Fehlerfunktion an den Stellen x i {\displaystyle x_{i}} alterniert ( sgn ⁡ ( f − p ) ( x i − 1 ) = − sgn ⁡ ( f − p ) ( x i ) {\displaystyle \operatorname {sgn}(f-p)(x_{i-1})=-\operatorname {sgn}(f-p)(x_{i})} ). Gibt es in zwei benachbarten Intervallen [ x i − 1 , x i ] {\displaystyle [x_{i-1},x_{i}]} und [ x i , x i + 1 ] {\displaystyle [x_{i},x_{i+1}]} jeweils nur eine Nullstelle s i ∈ [ x i − 1 , x i ] {\displaystyle s_{i}\in [x_{i-1},x_{i}]} beziehungsweise s i + 1 ∈ [ x i , x i + 1 ] {\displaystyle s_{i+1}\in [x_{i},x_{i+1}]} , dann ist die Fehlerfunktion im Intervall [ s i , s i + 1 ] {\displaystyle [s_{i},s_{i+1}]} nur nicht-negativ oder nur nicht-positiv. (Aber auch wenn es mehrere Nullstellen gibt, gilt das Weitere – die Extrema sind ggf. nur nicht so ausgeprägt.) Wegen der Stetigkeit der Fehlerfunktion gibt es nach dem Satz vom Minimum und Maximum in jedem Teilintervall [ s i , s i + 1 ] {\displaystyle [s_{i},s_{i+1}]} (lokale) Extremstellen x ~ i {\displaystyle {\tilde {x}}_{i}} . Im Ergebnis ist x ~ 0 := a {\displaystyle {\tilde {x}}_{0}:=a} das erste Extremum, die nachfolgenden Extrema alternieren zwischen Maximum und Minimum bis zum letzten Extremum x ~ n + 1 := b {\displaystyle {\tilde {x}}_{n+1}:=b} .

Nebenbei fällt die (absolute) Güte der Approximation in Form von max i | f ( x ~ i ) − p ( x ~ i ) | {\displaystyle \max _{i}|f({\tilde {x}}_{i})-p({\tilde {x}}_{i})|} an. Ist sie unbefriedigend, kann man mit der neu gewonnenen Referenz x ~ 0 , x ~ 1 , … , x ~ n + 1 {\displaystyle {\tilde {x}}_{0},{\tilde {x}}_{1},\ldots ,{\tilde {x}}_{n+1}} iterieren. Es kann aber auch interessant sein, den Grad n {\displaystyle n} der Approximation durch Einfügen von Zwischenpunkten in die Referenz zu erhöhen. Das Beispiel 2 im Artikel Alternantensatz zeigt, wie die Qualität der dortigen besten Approximation vom Polynomgrad abhängt.

Unter geeigneten Voraussetzungen, die Funktion f {\displaystyle f} betreffend, konvergiert das Verfahren lokal quadratisch.^[4]

• Leçons sur l’approximation des fonctions d’une variable réelle. Gauthier-Villars, Paris 1919, 1952; archive.org
• C. de Boor, A. Pinkus: Proof of the conjectures of Bernstein and Erdös concerning the optimal nodes for polynomial interpolation. In: Journal of Approximation Theory. 24. Jahrgang, 1978, S. 289–303, doi:10.1016/0021-9045(78)90014-X. 

1. ↑ E. Ya. Remez: Sur la détermination des polynômes d’approximation de degré donnée. In: Comm. Soc. Math. Kharkov, 1934, 10, S. 41.
   Sur un procédé convergent d’approximations successives pour déterminer les polynômes d’approximation. In: Compt. Rend. Acad. Sc., 1934, 198, S. 2063
   Sur le calcul effectiv des polynômes d’approximation de Tschebyscheff. In: Compt. Rend. Acad. Sc., 1934, 199, S. 337–340.
2. ↑ René Grothmann: Skriptum Approximationstheorie. (PDF) 1.69. Definition
3. ↑ Die angegebene Matrix hat Vandermonde-Matrizen als Untermatrizen.
   Die Lösung des Gleichungssystems kostet bei vielen Standardpaketen O ( n 3 ) {\displaystyle {\mathcal {O}}(n^{3})} , kann aber auch in O ( n 2 ) {\displaystyle {\mathcal {O}}(n^{2})} geschafft werden (s. Polynominterpolation#Newtonscher Algorithmus).
4. ↑ René Grothmann: Skriptum Approximationstheorie. (PDF) 1.71. Bemerkung