ラグランジュ補間

数値解析におけるラグランジュ補間（ラグランジュほかん、英: Lagrange interpolation）は、多項式補間に用いられる。相異なる点の集合 x_j および数値 y_j に対し、そのラグランジュ補間多項式は、各 x_j において対応する値として y_j をとるような次数最小の多項式である。このように次数最小の多項式は一意に決まるが、決定する方法は複数存在するため、「ラグランジュ補間多項式」という名称をその一意な多項式の「ラグランジュ形」というふうに言及するのは正確でない。

四点((−9, 5), (−4, 2), (−1, −2), (7, 9)) に対する三次補間多項式 L(x) (黒破線) は基底多項式 y₀ℓ₀(x), y₁ℓ₁(x), y₂ℓ₂(x), y₃ℓ₃(x) のスケールを変えたものの和である。この補間多項式は与えられた四つの制御点を通り、各スケールされた基底多項式は対応する制御点を通り、それ以外の三つの制御点に対応する x のところでは 0 になっている。

名称はジョゼフ＝ルイ・ラグランジュに因んだものだが、ラグランジュの発表する1795年よりも以前に、この方法を初めて発見したのは1779年の(エドワード・ワーリング)である。ラグランジュの結果はレオンハルト・オイラーが1783年に発表したより複雑な形の公式の簡単な帰結となるものであった^[1]

ラグランジュ補間多項式は数値積分法の一種ニュートン–コーツ法でも用いられ、また有限体上で計算されたラグランジュ補間多項式は暗号理論における(シャミアの秘密分散法)（英語版）でも用いられる。

ラグランジュ補間は巨大振幅に関するルンゲ現象の影響を受けやすい。また評価点 x_j の変更に関して補間の計算を全くやり直す必要があるから、そのような目的では変更が容易にできるニュートン補間がしばしば用いられる。

定義

k + 1 個のデータ点集合 ${\textstyle (x_{0},y_{0}),\ldots ,(x_{j},y_{j}),\ldots ,(x_{k},y_{k})}$  はどの二つの x_j も同じでないとする。これらのデータのラグランジュ型の補間多項式とは、ラグランジュ基底多項式

$${\displaystyle \ell _{j}(x):=\prod _{0\leq m\leq k \atop m\neq j}{\frac {x-x_{m}}{x_{j}-x_{m}}}={\frac {(x-x_{0})}{(x_{j}-x_{0})}}\cdots {\frac {(x-x_{j-1})}{(x_{j}-x_{j-1})}}{\frac {(x-x_{j+1})}{(x_{j}-x_{j+1})}}\cdots {\frac {(x-x_{k})}{(x_{j}-x_{k})}}}$$

 

の線型結合

$${\displaystyle L(x):=\sum _{j=0}^{k}y_{j}\ell _{j}(x)}$$

 

を言う。

まず最初の仮定により ${\textstyle x_{j}-x_{m}\neq 0}$  であるから、この式は常に矛盾なく定まる。${\textstyle x_{i}=x_{j}}$  かつ ${\displaystyle y_{i}\neq y_{j}}$  となるような対を許さない理由は、必要な補間多項式 L (${\textstyle y_{i}=L(x_{i})}$ ) が存在しないからである（各引数 x_i に対する値は一つでなければならない）。他方、$y_{i}=y_{j}$  もそれら二点が実際には単一の点となる。

任意の i ≠ j に対して、基底多項式 ℓ_j(x) は分子に (x − x_i) を因子として持つから、x = x_i のとき積全体も零となる。他方、i = j のときは明らかに ℓ_i(x_i) = 1 が成り立つ。すなわち、 ${\textstyle \ell _{j}(x_{i})=\delta _{ij}}$  である。ここに右辺はクロネッカーのデルタ。したがって各点 x_i において ${\textstyle L(x_{i})=y_{i}+0+0+\dotsb +0=y_{i}}$  となり L は所期の函数の補間を与えるものとなる。

注意

 

ラグランジュ多項式の集合に対する補間の発散の例

ラグランジュ型の補間では、多項式補間の線型性と補間多項式の一意性があるため、証明や理論的な議論では都合がよい。この一意性はヴァンデルモンド行列の可逆性（同じことだがヴァンデルモンド行列式が至る所消えないこと）からも導けることである。

しかしながら、構成からわかる通り、節点 x_k を変更するごとに、ラグランジュ基底多項式をすべて計算し直さなければならない。そこで実用上（あるいは計算上）の目的では、重心形のラグランジュ補間（後述）やニュートン補間を用いるほうが良い。

等間隔に取った節点に関するラグランジュ補間やその他の補間は、真の函数の上下に振動する多項式を生じるものである。この振る舞いは節点の数を多くするにつれてルンゲ現象と呼ばれる発散に逢着しやすくなっていく。この問題は、補間に用いる点を(チェビシェフ節点)（英語版）に選ぶことで解消できる^[2]。

線型代数学的説明

補間問題を解くことは、逆行列を計算する線型代数学の問題につながる。標準的な(単項式基底)を用いた補間多項式 ${\textstyle L(x)=\sum \limits _{j=0}^{k}x^{j}m_{j}}$  では、L(x_i) = y_i を L の係数 m_j に関してヴァンデルモンド行列 ${\textstyle ((x_{i})^{j})}$  を逆に解かなければならない。対して、ラグランジュ基底を用いて ${\textstyle L(x)=\sum \limits _{j=0}^{k}l_{j}(x)y_{j}}$  を作れば、この場合先ほどはヴァンデルモンド行列が現れた部分には単位行列 ${\textstyle (\delta _{ij})}$  が現れ、逆行列は単位行列自身であるから、ラグランジュ基底は自動的に逆に解かれていることになる。

この構成は中国の剰余定理と類似対応している（つまり、素数を法とする整数の剰余を調べる代わりに、一次式を法とする多項式の剰余を見ている）。

重心形

ラグランジュ基底多項式を

$${\displaystyle {\begin{aligned}\ell (x)&=(x-x_{0})(x-x_{1})\cdots (x-x_{k})\\[5pt]\ell '(x_{j})&={\frac {d\ell (x)}{dx}}{\Big |}_{x=x_{j}}=\prod _{i=0,i\neq j}^{k}(x_{j}-x_{i})\end{aligned}}}$$

 

を用いて、${\textstyle \ell _{j}(x)={\frac {\ell (x)}{\ell '(x_{j})(x-x_{j})}}}$ と書き直す。これは重心重み付け (barycentric weight)^[3]を ${\textstyle w_{j}={\frac {1}{\ell '(x_{j})}}}$  と定めれば簡潔に

$${\displaystyle \ell _{j}(x)=\ell (x){\frac {w_{j}}{x-x_{j}}}}$$

 

と書くことができる、これを重心ラグランジュ補間の「第一形」と呼ぶ。

この形の多項式補間を考えることの利点は、補間多項式 ${\textstyle L(x)=\ell (x)\sum \limits _{j=0}^{k}{\frac {w_{j}}{x-x_{j}}}y_{j}}$  を評価するときに、重み w_j が事前に分かっていれば、O(n) で計算できることである（ラグランジュ基底 ℓ_j(x) を個別に計算するのは O(n²) 掛かる）。もうひとつ重心形補間の利点として、新しい節点 x_k+1 の追加も各 w_j を ${\textstyle (x_{j}-x_{k+1})}$  で割って、新しい w_j+1 を計算するだけで容易にできる点が挙げられる。

さらに第一形を単純化することもできて、まず定数函数 g(x) ≡ 1 の重心補間 ${\textstyle g(x)=\ell (x)\sum \limits _{j=0}^{k}{\frac {w_{j}}{x-x_{j}}}}$  を計算してから L を g で割れば、得られる

${\displaystyle L(x)={\frac {\textstyle \sum \limits _{j=0}^{k}{\frac {w_{j}}{x-x_{j}}}y_{j}}{\sum \limits _{j=0}^{k}{\frac {w_{j}}{x-x_{j}}}}}}$  は与えられた節点における補間性を失わない。この補間多項式を重心補間多項式の「第二形」あるいは「真の形」という。真の重心補間多項式は、L の各節点における評価に際してラグランジュ基底 ℓ を評価しなくてよいという点で有利である。

誤差項

与えられた函数 f を、節点 x₀, …, x_n において次数 n の多項式で補完するとき、誤差 ${\displaystyle R(x)=f(x)-L(x)}$  は

$${\displaystyle R(x)=f[x_{0},\ldots ,x_{n},x]\ell (x)=\ell (x){\frac {f^{(n+1)}(\xi )}{(n+1)!}},\quad \quad x_{0}<\xi <x_{n},}$$

 

と表せる^[4]。ただし、${\textstyle f[x_{0},\ldots ,x_{n},x]}$  は差商である^{[注釈 1]}. またこの誤差項を複素領域における周回積分

$${\displaystyle R(z)={\frac {\ell (z)}{2\pi i}}\int _{C}{\frac {f(t)}{(t-z)(t-z_{0})\cdots (t-z_{n})}}dt={\frac {\ell (z)}{2\pi i}}\int _{C}{\frac {f(t)}{(t-z)\ell (t)}}dt}$$

 

と書くこともできる。

この誤差項によって、誤差の範囲を

$${\displaystyle |R(x)|\leq {\frac {(x_{n}-x_{0})^{n+1}}{(n+1)!}}\max _{x_{0}\leq \xi \leq x_{n}}|f^{(n+1)}(\xi )|}$$

 

と見積もることができる。

脚注

1. ^ 差分商とも呼ぶ。近い概念として、有限差分の商 Δy/Δx を差分商あるいは差商と呼ぶ場合もあるので混同してはならない。

参考文献

1. ^ Meijering, Erik (2002), “A chronology of interpolation: from ancient astronomy to modern signal and image processing”, Proceedings of the IEEE 90 (3): 319–342, doi:10.1109/5.993400 .
2. ^ Quarteroni, Alfio; Saleri, Fausto (2003), Scientific Computing with MATLAB, Texts in computational science and engineering, 2, Springer, p. 66, ISBN (9783540443636), https://books.google.com/books?id=fE1W5jsU4zoC&pg=PA66 .
3. ^ Jean-Paul Berrut & (Lloyd N. Trefethen) (2004). “Barycentric Lagrange Interpolation”. (SIAM) Review 46 (3): 501-517. doi:10.1137/S0036144502417715. 
4. ^ Abramowitz and Stegun, "Handbook of Mathematical Functions," p.878

• 高橋大輔 『数値計算』 岩波書店、1996年。

関連項目

• ネヴィルのアルゴリズム
• ニュートン補間: ニュートン形の補間多項式
• バーンスタイン多項式: バーンスタイン形の補間多項式
• (カールソンの定理)（英語版）
• (ルベーグ定数 (補間法))（英語版）
• (Chebfun)（英語版）
• (ニュートン級数の一覧)（英語版）
• (フロベニウス共変行列)（英語版）
• (シルベスターの公式)（英語版）: ラグランジュ–シルベスター補間

外部リンク

• 『ラグランジュの補間公式とその応用例』 - 高校数学の美しい物語
• Hazewinkel, Michiel, ed. (2001), "Lagrange interpolation formula", Encyclopaedia of Mathematics, Springer, ISBN 978-1-55608-010-4。
• Weisstein, Eric W. "Lagrange Interpolating Polynomial". MathWorld (英語).
• (Lagrange Interpolation Formula) at ProofWiki
• (Lagrange Polynomial Approximation) at ProofWiki

• ALGLIB has an implementations in C++ / C# / VBA / Pascal.
• GSL has a polynomial interpolation code in C
• SO has a MATLAB example that demonstrates the algorithm and recreates the first image in this article
• Lagrange Method of Interpolation — Notes, PPT, Mathcad, Mathematica, MATLAB, Maple at Holistic Numerical Methods Institute
• Lagrange interpolation polynomial on www.math-linux.com
• Dynamic Lagrange interpolation with JSXGraph
• Numerical computing with functions: The Chebfun Project
• Worksheet Function for Bicubic Lagrange Interpolation
• Lagrange polynomials in Python