多項式補間のソースを表示

'''多項式補間'''（たこうしきほかん、{{lang-en-short|polynomial interpolation}}）は、[[数値解析]]において、与えられたデータ群を[[多項式]]で[[内挿]]（補間）することである。言い換えれば、[[標本調査]]などで得たデータ群について、それらを正確に通る多項式を見つけることである。

== 用途 ==
多項式をより複雑な曲線の近似として使う場合もあり、[[タイポグラフィ]]における文字の形状をいくつかの点で表すなどの例がある。関連する用途としては、[[自然対数]]や[[三角関数]]の値を求める際に、[[数表]]に掲載されている点から多項式補間で必要な値を求める場合がある。これは特定の値を直接求めようとするよりもかなり高速に計算できる。多項式補間はまた、[[数値積分]]や[[常微分方程式]]の数値解を求めるアルゴリズムの基盤にもなっている。

多項式補間は、[[Karatsuba法]]や{{仮リンク|Toom-Cook法|en|Toom–Cook multiplication}}といった乗算アルゴリズムの基盤であり、積を定義する多項式上の点間の補間が積自体を生成する。例えば ''a'' = ''f''(''x'') = ''a''<sub>0</sub>''x''<sup>0</sup> + ''a''<sub>1</sub>''x''<sup>1</sup> + ... と ''b'' = ''g''(''x'') = ''b''<sub>0</sub>''x''<sup>0</sup> + ''b''<sub>1</sub>''x''<sup>1</sup> + ... があるとき、積 ''ab'' は ''W''(''x'') = ''f''(''x'')''g''(''x'') と等しい。''f''(''x'') と ''g''(''x'') が小さいときの ''x'' における ''W''(''x'') 上の点を求めることで、その曲線上の点が得られる。それらの点に基づく補間によって ''W''(''x'') の項とさらには積 ''ab'' が得られる。Karatsuba法の場合、穏当なサイズの入力であっても単純な乗算より高速である。特にハードウェアの並列性を利用するとその特性が発揮される。

== 定義 ==
''n''+1 個のデータ点 (''x''<sub>''i''</sub>,''y''<sub>''i''</sub>) があり、''x''<sub>''i''</sub> には同じ値のものがないとする。ここで、

:<math>p(x_i) = y_i,\;  i=0,\ldots,n</math>

という属性を持つ最大 ''n''次の多項式を求める。[[w:Unisolvent functions|unisolvence]]定理によれば、このような多項式 ''p'' は必ず存在し、一意に定まる。

より正確に言えば、その定理は ''n''+1 個の補間ノード (''x''<sub>''i''</sub>) について、多項式補間は次の線型[[全単射]]を定義する。
:<math>L_n:\mathbb{K}^{n+1} \to \Pi_n</math>
ここで <math>\Pi_n</math> は次数が ''n'' 以下の多項式の[[ベクトル空間]]である。

== 補間多項式の構築 ==
[[画像:Interpolation_example_polynomial.svg|thumb|right|赤い点はデータ点 (''x''<sub>''k''</sub>,''y''<sub>''k''</sub>)T を表し、青い曲線は補間多項式を表している。]]
補間多項式が次のような形式であるとする。
:<math>p(x) = a_n x^n + a_{n-1} x^{n-1} + \cdots + a_2 x^2 + a_1 x + a_0 \qquad (1) </math>
''p'' がデータ点群を補間するとは、次を意味する。
:<math>p(x_i) = y_i \qquad\mbox{for all } i \in \left\{ 0, 1, \dots, n\right\}</math>
ここで式 (1) に置換すると、係数 <math>a_k</math> の[[線型方程式系]]が得られる。これを行列とベクトルで表すと次のようになる。
:<math>\begin{bmatrix}
x_0^n & x_0^{n-1} & x_0^{n-2} & \ldots & x_0 & 1 \\
x_1^n & x_1^{n-1} & x_1^{n-2} & \ldots & x_1 & 1 \\
\vdots & \vdots & \vdots & & \vdots & \vdots \\
x_n^n & x_n^{n-1} & x_n^{n-2} & \ldots & x_n & 1 
\end{bmatrix}
\begin{bmatrix}
a_n \\
a_{n-1} \\
\vdots \\
a_0
\end{bmatrix}
=
\begin{bmatrix}
y_0 \\
y_1 \\
\vdots \\
y_n
\end{bmatrix}
</math>
補間多項式 <math>p(x)</math> を構築するには、この方程式系を <math>a_k</math> について解かなければならない。

左辺の行列を[[ファンデルモンド行列]]と呼ぶ。その[[行列式]]はゼロではなく、唯一の補間多項式が存在する。

ファンデルモンド行列の[[条件数]]は大きいため<ref>{{cite journal|last=Gautschi|first=Walter|title=Norm Estimates for Inverses of Vandermonde Matrices|journal=Numerische Mathematik|volume=23|pages=337–347|year=1975|doi=10.1007/BF01438260}}</ref>、係数 <math>a_i</math> を求めるのに[[ガウスの消去法]]を使うと、大きな誤差を生じる。そのため、<math>\mathcal O(n^3)</math> を要する[[ガウスの消去法]]の代わりに、ファンデルモンド行列の特性を利用した <math>\mathcal O(n^2)</math> の[[数値的安定性|数値的に安定な]]解法がいくつも提案されている<ref>{{cite journal|last=Higham|first=N. J.|title=Fast Solution of Vandermonde-Like Systems Involving Orthogonal Polynomials|journal=IMA Journal of Numerical Analysis|volume=8|pages=473–486|year=1988|doi=10.1093/imanum/8.4.473}}</ref><ref>{{cite journal|last=Björck|first=Å|coauthors=V. Pereyra|title=Solution of Vandermonde Systems of Equations|journal=Mathematics of Computation|volume=24|number=112|pages=893–903|year=1970|doi=10.2307/2004623}}</ref><ref>{{cite journal|author=Calvetti, D and Reichel, L|title=Fast Inversion of Vanderomnde-Like Matrices Involving Orthogonal Polynomials|journal=BIT|number=33|pages=473–484|year=1993|doi=10.1007/BF01990529|volume=33}}</ref>。これらの手法はまず多項式の[[ニュートン補間]]を構築し、それを上掲のような単項形式に変換する。

== ファンデルモンド以外の解法 ==
ここで構築しようとしているのは、''n''次の多項式空間であるベクトル空間 <math>\Pi_n</math> にある唯一の補間多項式である。<math>\Pi_n</math> の{{仮リンク|単項式基底|en|Monomial basis}}を用いると、ファンデルモンド行列を解いて補間多項式の係数群 <math>a_k</math> を構築しなければならない。これは（コンピュータのクロック数で数えると）非常にコストがかかる操作である。<math>\Pi_n</math> の別の基底を選択すると係数の計算を単純化できるが、補間多項式を単項式基底で表現するために追加の計算が必要になる。

1つの手法として、補間多項式を[[ニュートン補間]]形式で書き、[[差商|差分商]]の手法で係数を構築する（これを{{仮リンク|ネヴィル法|en|Neville's algorithm}}という）。そのコストは [[ランダウの記号|O]]<math>(n^2)</math> で、一方ガウスの消去法のコストは O<math>(n^3)</math> である。さらに、点が追加されたときそれまでの計算結果に若干加える形でよいという利点がある（他の手法では、全部の計算をやり直す必要がある）。

別の手法として、[[ラグランジュ補間]]形式の補間多項式を使う方法がある。得られる方程式から上の定義に指定されている条件に当てはまる補間多項式を即座に示すことができる。
多項式補間を行なうのに便利な方法として、重心形式補間法(barycentric interpolation)がある<ref>Jean-Paul Berrut and Lloyd N. Trefethen: "Barycentric Lagrange Interpolation", SIAM Review, Vol.46, No.3, pp.501-517 (2004).</ref>。これは<math>n</math>個の分点があらかじめ指定されているとき，補間1回あたり計算量<math>O(n)</math>で関数の補間多項式の値が求まる。

[[バーンスタイン多項式|ベルンシュテイン多項式]]（{{en|Bernstein polynomials}}）は[[セルゲイ・ベルンシュテイン]]が[[ストーン＝ワイエルシュトラスの定理|ワイエルシュトラスの近似定理]]の構築的証明に使ったが、現在ではコンピュータグラフィックスでよく使われる[[ベジェ曲線]]の元になっている。(注：ベルンシュテイン多項式は、端点以外では分点上で関数値は一致しないので、「多項式補間」ではないが、次数を上げるとともに元の関数に区間内で一様に収束する多項式の列を与える近似法である。)

== 補間誤差 ==
関数 ''f'' を ''x''<sub>''0''</sub>,...,''x''<sub>''n''</sub> というノードを結ぶ''n''次多項式で補間したとき、その誤差は以下のようになる。

:<math>f(x) - p_n(x) = f[x_0,\ldots,x_n,x] \prod_{i=0}^n (x-x_i) </math>

ここで
:<math>f[x_0,\ldots,x_n,x]</math>
は差分商を表す。''f'' がノード ''x''<sub>''i''</sub> と ''x'' を含む微小な区間 ''I'' において ''n''+1 階連続微分可能ならば、''I'' 上のある <math>\xi</math> について誤差をラグランジュ形式で次のように書くことができる。

:<math> f(x) - p_n(x) = \frac{f^{(n+1)}(\xi)}{(n+1)!} \prod_{i=0}^n (x-x_i) </math>

したがって、[[テイラーの定理]]のラグランジュ形式の剰余項は、全ての補間ノード ''x''<sub>''i''</sub> が同一であるような特殊ケースの補間誤差である。

<math>x_i = x_0 + ih</math> のように等間隔な補間ノードの場合、補間誤差は O<math>(h^n)</math> になる。しかし、このことから <math>n \to \infty</math> となったときに何が起きるかはわからない。それについては後述の[[#収束属性|「収束属性」の節]]で扱う。

以上から、積 | &prod; (''x'' &minus; ''x''<sub>''i''</sub>) | を可能な限り小さくするような補間点 ''x''<sub>''i''</sub> の選択がありうることがわかる。これを達成しているのが{{仮リンク|チェビシェフノード|en|Chebyshev nodes}}である。

== ルベーグ定数 ==
補間点 ''x''<sub>0</sub>, ..., ''x''<sub>''n''</sub> を設定し、全ての補間ノードが区間 [''a'', ''b''] にあるとする。補間とはすなわち、関数 ''f'' から多項式 ''p'' への写像をする過程である。これにより、[''a'', ''b''] 上の全ての連続関数の空間 ''C''([''a'', ''b'']) からそれ自身への写像 ''X'' を定義する。写像 ''X'' は線型であり、''n''次以下の多項式の部分空間 &Pi;<sub>''n''</sub> 上の[[射影]]である。

ルベーグ定数 ''L'' は ''X'' の[[ノルム#無限次元ベクトル空間のノルム|作用素ノルム]]として定義される。それは次のようになる（[[:en:Lebesgue's lemma|ルベーグの補題]]の特殊ケース）。
:<math> \|f-X(f)\| \le (L+1) \|f-p^*\| </math>
言い換えれば、補間多項式は最良の近似に比べて最大で係数 (''L''+1) のぶんだけ悪い。従って補間ノード集合は ''L'' が小さいほどよいことを示唆している。特にチェビシェフノードは次のようになる。
:<math> L \ge \frac2\pi \log(n+1) + C \quad\mbox{for some constant } C. </math>
つまり、チェビシェフノードは多項式補間には最適の選択であり、同程度の良さを等間隔ノードで実現しようとすると ''n'' を指数関数的に成長させる必要がある。しかし、それらノードは最適ではない。

== 収束属性 ==
補間多項式の次数 ''n'' を無限大に漸近させたとき、その一連の補間多項式が補間対象の関数に収束するような関数の種類や補間ノードのシーケンスはどのようなものだろうか? 収束は他にも、様々なノルムについて考えられる。

等間隔ノードでは状況は悪く、無限に微分可能な関数では一様収束も保証されない。[[ルンゲ現象|カール・ルンゲの古典的例]]として、区間 [&minus;5, 5] の関数 ''f''(''x'') = 1 / (1 + ''x''<sup>2</sup>) がある。この場合、補間誤差 ||''f'' &minus; ''p''<sub>''n''</sub>||<sub>&infin;</sub> は ''n'' &rarr; &infin; に従って大きくなる。もう1つの例として、区間 [&minus;1, 1] の関数 ''f''(''x'') = |''x''| がある。この場合の補間多項式は ''x'' = &minus;1, 0, 1 の3点以外では各点収束しない<ref>Watson (1980, p. 21) attributes the last example to Bernstein (1912).</ref>。

補間ノードの選択によって収束属性がよくなるだろうか? 次の'''定理'''がある程度答になるだろう。

: 任意の関数 ''f''(''x'') が区間 [''a'',''b''] 上で連続であるとき、一連の補間多項式 <math>p_n(x)</math> が [''a'',''b''] 上で ''f''(''x'') に一様に収束するようなノードの並びが存在する。

証明は次の通り。最良の近似の一連の多項式 <math>p^*_n(x)</math> が ''f''(''x'') に一様に収束することは明らかである（ワイエルシュトラスの近似定理より）。したがって、ここでは個々の <math>p^*_n(x)</math> をあるノード群の補間によって得られることを示せばよい。しかし、チェビシェフの交替定理として知られる最良近似の多項式ではそれが真である。特に補間多項式は ''f''(''x'') と少なくとも ''n''+1 回交差する。交差点が補間ノードとなるよう選択すると、最良の近似多項式と同一の補間多項式が得られる。

しかしこの方法の欠陥は、関数 ''f''(''x'') 毎に改めて補間ノードを計算しなければならず、しかもそのアルゴリズムは数値的に実装するのが難しい点である。どんな連続関数 ''f''(''x'') についても一連の補間多項式が収束するようなノードの組合せがあるだろうか? 答は残念ながら次の'''定理'''で示される通り否定的である。

: ノードの任意の組合せについて、区間 [''a'',''b''] において一連の補間多項式が発散するような連続関数 ''f''(''x'') が存在する<ref>Watson (1980, p. 21) attributes this theorem to Faber (1914).</ref>。

その証明はルベーグ定数の下限推定を使うもので、それは上で定義したように ''X''<sub>''n''</sub>（すなわち &Pi;<sub>''n''</sub> 上の射影作用素）の作用素ノルムである。ここで、次のような性質のノードの組合せを探す。
:<math>\lim_{n \to \infty} X_n f = f,</math> for any  <math>f \in C([a,b]) </math>
{{仮リンク|一様有界性原理|en|Uniform boundedness principle}}によると、''X''<sub>''n''</sub> のノルムが一様有界であるときだけこのようなノード群が存在するが、<math>\|X_n\|\geq \frac{2}{\pi} \log(n+1)+C</math> であることから、そのようなノード群は存在しない。

例えば、等間隔点を補間ノードとして選択したとき、[[ルンゲ現象]]による発散が見られる。なお、その関数は単に連続であるだけでなく、区間 [&minus;1, 1] において無限に微分可能である。チェビシェフノードより良い例は、次の'''定理'''により、見つけるのがさらに難しい。

: 区間 [&minus;1, 1] で[[絶対連続]]なあらゆる関数について、チェビシェフノードで構築した一連の補間多項式は一様に ''f''(''x'') に収束する。

== 関連する概念 ==
[[ルンゲ現象]]とは、''n'' が大きくなると補間多項式がデータ点間でより大きく振動するようになる現象である。この問題は一般に[[スプライン補間]]によって解決する。そのときの補間曲線は多項式ではなく[[スプライン曲線]]であり、低次のいくつかの多項式の連鎖になっている。

[[調和解析|調和関数]]で[[周期関数]]を補間する場合、一般に[[フーリエ級数]]を用い、例えば[[離散フーリエ変換]]などで行っている。これは調和基底関数による多項式補間と見ることもできる。

{{仮リンク|エルミート補間|en|Hermite interpolation}}問題は、多項式 ''p'' の値が与えられるだけでなく、いくつかの導関数値も与えられる場合である。[[:en:Birkhoff interpolation|バーコフ補間]]は、それをさらに一般化し、''p'' の値そのものは与えられず微分値だけを与えられる場合も含む。

微分方程式や積分方程式を解く{{仮リンク|コロケーション法|en|Collocation_method}}は、多項式補間に基づいている。

[[:en:Polynomial and rational function modeling|有理関数モデリング]]の技法は、多項式関数の比を考慮した一般化である。

== 脚注 ==
{{Reflist}}

== 参考文献 ==
* Kendell A. Atkinson (1988). ''An Introduction to Numerical Analysis'' (2nd ed.), Chapter 3. John Wiley and Sons. ISBN 0-471-50023-2.
* Sergei N. Bernstein (1912), Sur l'ordre de la meilleure approximation des fonctions continues par les polynômes de degré donné. ''Mem. Acad. Roy. Belg.'' '''4''', 1&ndash;104.
* L. Brutman (1997), Lebesgue functions for polynomial interpolation — a survey, ''Ann. Numer. Math.'' '''4''', 111&ndash;127.
* Georg Faber (1912), Über die interpolatorische Darstellung stetiger Funktionen, ''Deutsche Math. Jahr.'' '''23''', 192&ndash;210.
* M.J.D. Powell (1981). ''Approximation Theory and Methods,'' Chapter 4. [[ケンブリッジ大学出版局|Cambridge University Press]]. ISBN 0-521-29514-9.
* Michelle Schatzman (2002). ''Numerical Analysis: A Mathematical Introduction,'' Chapter 4. [[:en:Clarendon Press|Clarendon Press]], Oxford. ISBN 0-19-850279-6.
* Endre S&uuml;li and David Mayers (2003). ''An Introduction to Numerical Analysis,'' Chapter 6. [[:en:Cambridge University Press|Cambridge University Press]]. ISBN 0-521-00794-1.
* G. Alistair Watson (1980). ''Approximation Theory and Numerical Methods''. John Wiley. ISBN 0-471-27706-1.

== 外部リンク ==
* [http://demonstrations.wolfram.com/InterpolatingPolynomial/ Interpolating Polynomial] by [[スティーブン・ウルフラム|Stephen Wolfram]], [[ウルフラム・リサーチ|The Wolfram Demonstrations Project]].

{{DEFAULTSORT:たこうしきほかん}}

[[Category:補間]]
[[Category:数学に関する記事]]